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The transition from stick to slip at a dry frictional interface occurs through the breaking of 
the junctions between the two contacting surfaces. Typically, interactions between the junctions 
through the bulk lead to rupture fronts propagating from weak and/or highly stressed regions, whose 
junctions break first. Experiments find rupture fronts ranging from quasi-static fronts with speeds 
proportional to external loading rates, via fronts much slower than the Rayleigh wave speed, and 
fronts that propagate near the Rayleigh wave speed, to fronts that travel faster than the shear wave 
speed. The mechanisms behind and selection between these fronts are still imperfectly understood. 
Here we perform simulations in an elastic 2D spring-block model where the frictional interaction 
between each interfacial block and the substrate arises from a set of junctions modeled explicitly. 
We find that a proportionality between material slip speed and rupture front speed, previously 
reported for slow fronts, actually holds across the full range of front speeds we observe. We revisit a 
mechanism for slow slip in the model and demonstrate that fast slip and fast fronts have a different, 
inertial origin. We highlight the long transients in front speed even in homogeneous interfaces, and 
we study how both the local shear to normal stress ratio and the local strength are involved in the 
selection of front type and front speed. Lastly, we introduce an experimentally accessible integrated 
measure of block slip history, the Gini coefficient, and demonstrate that in the model it is a good 
predictor of the history-dependent local static friction coefficient of the interface. These results will 
contribute both to building a physically-based classification of the various types of fronts and to 
identifying the important mechanisms involved in the selection of their propagation speed. 

PACS numbers: 81.40.Pq, 46.55.+d, 62.20.Qp, 46.50.+a 


I. INTRODUCTION 

The onset of sliding at a frictional interface occurs 
through the breaking of the many contacts that were 
preventing the relative motion of the surfaces. When 
a single contact breaks, the stress it bore is redistributed 
to its neighbors, bringing them closer to or past their 
load-bearing capacity. In extended frictional interfaces 
(i.e. larger than a characteristic correlation length scale), 
this process can lead to propagating ruptures - rupture 
fronts. The recent direct observation of rupture fronts 
in laboratory friction experiments (see e.g. [1-5]) have 
provided new insights and opened new questions. It was 
found, for instance, that not all fronts span the entire in¬ 
terface [6, 7]. The selection of the propagation length of 
the fronts has been intensely investigated [8-14]. It was 
also found that the fronts can propagate at a variety of 
speeds, either quasi-statically [4, 5, 15], at speeds close 
to that of surface waves (sub-Rayleigh) [3, 16, 17], at 
speeds faster than the shear wave speed Cg (supershear) 
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[1, 3, 16-18], or at speeds one or two orders of magnitude 
slower than the Rayleigh speed (slow) [3, 16, 19, 20]. In 
this context, the present paper is mainly devoted to the 
study of the mechanisms responsible for front speed se¬ 
lection. 

Experiments to shed light on the nature of the rup¬ 
ture fronts have been performed. The authors of [16, 21] 
placed arrays of sensors close to the interface and used 
continuum theory to infer the properties of the elastic 
fields at the interface. Svetlizky and Fineberg [22] re¬ 
cently showed that the dynamic fields associated with 
sub-Rayleigh fronts are consistent with the ones pre¬ 
dicted by linear elastic fracture mechanics [23]. The au¬ 
thors of [5, 20, 24] used microtextured surfaces where 
each contact can be tracked individually to follow the 
rupture fronts directly at the individual micro-contact 
level. Common to all the experiments mentioned till 
now is that they seek to increase both the spatial and 
the temporal resolution of the behavior of the very inter¬ 
face. These experiments can be usefully complemented 
by computer simulations, which by their nature pro¬ 
vide complete information of all quantities in the mod¬ 
els they implement, thus filling in the information gaps 
that remain despite the experimental progress; for ex¬ 
ample, simulations can provide simultaneous access to 
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shear stress, normal stress, local contact strength and 
front propagation, a combination which remains hard to 
access experimentally. 

To represent the propagation of the front separating a 
stuck part of the interface from a slipping one, models of 
the transition to sliding need to include at least one level 
of discretization of the macroscopic interface. Depending 
on the model, this so-called mesoscopic scale may itself 
host a population of smaller scaled micro-junctions, such 
as the micro-contacts forming the multi-contact between 
rough surfaces (see e.g. [25-27]). However, the emerg¬ 
ing collective behavior of the micro-junctions can also be 
lumped into a local friction law, often inspired by laws de¬ 
veloped empirically from macroscale experiments, acting 
directly at the meso-scale (see e.g. [7-9, 28-35]). The 
existing models also differ in how they treat the stress 
transfer in the bulk. The options range from ID sys¬ 
tems [7, 8, 10, 12, 25, 26, 28-30, 34, 35] via 2D models 
[9, 11, 27, 31-33] to a full 3D visco-elasto-plastic dis¬ 
cretization. In principle, the most comprehensive fric¬ 
tion models could be combined with the most compre¬ 
hensive bulk models, but the dynamics that arise from 
this approach tend to be nearly as complicated and dif¬ 
ficult to disentangle as in the experiments themselves, 
and instead, authors have focused on one or a few model 
properties. 

Continuum models of friction, from the Amontons- 
Coulomb description to more sophisticated rate-and- 
state friction laws [36-40], are successful when they re¬ 
produce a robust average behavior of the myriad micro¬ 
junctions that make up each mesoscopic region of the fric¬ 
tional interface. However, by their nature, these models 
do not explain how the individual micro-junctions evolve 
and interact to produce the overall friction behavior. To 
approach this question, numerical [41-43] and analytical 
[44-47] models have been made that explicitly include a 
set of junctions, each representing one or a few micro¬ 
scopic contacts. The main missing ingredient in these 
models is a solid foundation for the junction evolution 
laws. In principle, this could be addressed by molecular 
dynamics simulations [48], but these simulations enable 
too short time scales and too small length scales for a sys¬ 
tematic study of the onset of macroscopic sliding. Even 
models in which junctions represent thermally activated 
single molecular bonds in a simplified way [49, 50] would 
require a very large number of junctions and prohibitive 
calculation times in order to resolve the spatiotemporal 
dynamics of the transition to sliding observed experimen¬ 
tally. 

The simplest bulk model that includes spatial and tem¬ 
poral structure in the transition to sliding is probably the 
ID spring-block model. Friction has been of interest to 
the earthquake community at least since Brace and By- 
erlee [51] suggested that shallow earthquakes are the slip 
events in a stick-slip cycle, and since the spring-block 
experiment and simulation by Burridge and Knopoff [52] 
the ID spring-block model has been popular in the earth¬ 
quake literature, see e.g. [53, 54]. As spatiotemporal data 


for the onset of sliding became available in laboratory 
friction experiments, the ID spring-block model was also 
applied in the friction literature [7, 10, 12, 25, 26, 30, 55]. 
The main limitations of the ID spring-block model are 
its inability to accurately reproduce the stress fields that 
arise in the experiments and the lack of a physical length 
scale [56, 57] unless such a length scale is introduced in 
the friction law (e.g. [10]). 

To better reproduce the experimental loading condi¬ 
tions and how they translate to heterogeneous shear and 
normal stress fields at the frictional interface, we used, 
in [9], a 2D spring-block model, which can be shown 
to reproduce 2D linear elasticity [58]. With Amontons- 
Coulomb friction this model agreed well with experi¬ 
ments for static measures related to the onset of sliding, 
such as for the length of precursors and the evolution 
of the normal stress along the interface, but the model 
did not capture the full dynamics of the rupture fronts. 
Radiguet et al. [11] studied the memory of the stress state 
through the passage of multiple ruptures in a visco-elastic 
2D finite element model with a slip-weakening friction 
law. They too focused on successive precursors rather 
than the dynamics of each rupture event. Kammer et al. 

[32] studied the properties of fast rupture front propaga¬ 
tion in a 2D finite element model with a static+velocity 
weakening dynamic friction law. Otsuki and Matsukawa 

[33] studied how the normal force and the size of the in¬ 
terface influence the effective macroscopic static friction 
coefficient, using a 2D finite element model with a lo¬ 
cal velocity-weakening Amontons-Coulomb friction law. 
While the above-mentioned models were able to model 
various aspects of the properties of rupture fronts, they 
do not provide a framework that within the same model 
is able to account for the full richness of front dynam¬ 
ics that was observed experimentally. In particular, the 
self-selection of front type leading to sub-Rayleigh, slow 
and super shear fronts, as well as the transitions between 
them, was missing. 

In [27] we combined an asperity model of the friction 
at the interface with a 2D elastic solver in order to ac¬ 
curately reproduce the experimental loading conditions 
used in [3]. In the asperity model we included a time 
scale inspired by the time scale identified in the same 
experimental system [59] that in the model controls the 
healing of the interface back to a fully pinned state after 
slipping. This combined model produced spatio-temporal 
features of the rupture dynamics very similar to those ob¬ 
served in the experiment. In particular, we reproduced 
the abrupt transitions between fast and slow front propa¬ 
gation, which can be understood from the underlying fast 
and slow slip mechanisms. Here, we use this very same 
model to gain significant additional insights into the re¬ 
lation between the micro-scale junction dynamics, the 
meso-scale slow and fast slip dynamics, and the macro¬ 
scale friction dynamics and front propagation. Our re¬ 
sults show that while fast slip and fast fronts are inertially 
controlled, slow slip and slow fronts are non-inertial and 
instead depend on properties of the friction model. Nev- 
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ertheless, the same proportionality relates front speed 
and slip speed in both velocity regimes. This suggests 
that the key to understanding complicated experimen¬ 
tal front dynamics lies in understanding the underlying 
slip dynamics. Further, we show that the spatial extent 
of transients in front speed is comparable to the size of 
our system (we take the system size from experiments). 
This suggests that, although front speed is influenced 
strongly by local quantities like the stress state and the 
local strength [27], front speed selection is intrinsically a 
non-local quantity which depends on the rupture history 
of the whole interface. Finally, we show that by applying 
the Gini coefficient, which is an integrated measure of in¬ 
equality commonly used in demography [60, 61, p. 186], 
we can characterize the complicated slip dynamics of a 
rupture event in the simulation and predict the subse¬ 
quent (slip-history-dependent) local frictional strength of 
the interface. 

The paper is structured as follows. In Section II we 
introduce the model. We first describe the motivating 
physical picture, and then make the model explicit by 
defining its equations and parameters. Then, we high- 
iight some simiiarities and differences with other popular 
friction models. Next follow four results sections. In Sec¬ 
tion III we show that the friction law introduced on the 
junction-level leads to stick-slip on the system level. We 
also revisit the slow slip mechanism identified in [27]. In 
Section IV we present simulations with both fast-only 
and fast-slow-fast rupture events. We identify a possi¬ 
ble signature of slow fronts in the macroscopic loading 
curve, and demonstrate the inertial nature of fast slip 
and fast fronts. We show how changes to parameters in 
the model or the state of the interface before rupture can 
turn fast-slow-fast fronts into fast-only fronts and vice 
versa. We map out the initial conditions that lead to 
fast-only and to slow fronts. In Section V we study the 
range of front speeds within each event and within each 
type of front (fast or slow). We show that in the model, 
the spatial extent of front speed transients is comparable 
to the current system size. We show that the influence 
of local stress and strength on front speed within each 
type of front is similar to their role in selecting the front 
type. We show that fast front speed is proportional to 
fast slip speed and that the constant of proportionality is 
the same as for slow fronts vs slow slip, providing a uni¬ 
versal relationship across the velocity scales in the model. 
In Section VI we come full circle by resolving how a rup¬ 
ture event sets up the state of the interface, which in 
turn determines the properties of the next rupture event. 
This links back to our results in [47]. Section VII is a 
brief discussion. The appendices provide additional de¬ 
tails we deemed important to the understanding of the 
model, but that would hamper the flow of the arguments 
were they to appear within the main text. These include 
the stencil used for finding the front speed numerically, 
a data collapse demonstrating the utility of some of our 
dimensionless parameters, and a detailed account of how 
we apply initial and boundary conditions. 


II. MODEL DESCRIPTION 

The frictional stability of a system made of two solids 
in contact locally depends on the level of normal and 
shear stresses at the contact interface. These interfa- 
cial stresses result from the external forces applied at the 
boundaries of the solids, transmitted through the bulk. 
Slip motion will in general be triggered when the local 
inter facial shear stress reaches a threshold, the level of 
which crucially depends on the interface behavior law at 
the microscale. The model we employ here is the same as 
in [27]. In this section we begin by describing the physi¬ 
cal aspects that underlie the model assumptions; we then 
describe the model and its parameters in detail; we end 
by discussing briefly the relationship to other models. 


A. Physical aspects 

The net contact between two solids generically con¬ 
sists of a large number of stress bearing micro-junctions. 
The properties of these junctions depend on the type of 
interface. For rough solids, each micro-junction corre¬ 
sponds to a micro-contact between asperities on the op¬ 
posing surfaces, whereas for smoother surfaces the junc¬ 
tions can be solidified patches of an adsorbate layer [43]. 
We include in our model the following three physical as¬ 
pects of the junction behavior. I) A micro-junction in 
its pinned state behaves elastically and can bear a shear 
force /t, provided /t remains smaller than a threshold 
/thres- When /thres IS reached, a local fracture-like event 
occurs, and the junction enters a slipping state. 2) In 
the slipping state, the micro-junction moves with the 
slider’s surface. The physical picture can be e.g. the 
micro-slipping of micro-asperities in contact or the flu¬ 
idization of an adsorbate layer. During slip, the micro¬ 
junction sustains some residual force /t = /slip, with /sUp 
smaller than /thres- 3) Slipping micro-junctions have a 
certain probability to disappear or relax. For example, a 
micro-contact disappears when an asperity moves away 
from its antagonist asperity by a typical distance equal to 
the mean size of micro-contacts, as classically considered 
for slow frictional sliding, e.g. in rate-and-state friction 
laws. However, another picture may arise from the sud¬ 
den release of energy when pinned junctions break. This 
energy will dissipate as heat in the region around the 
micro-j unction [59]. The rise in temperature will signif¬ 
icantly increase the rate of a thermally activated relax¬ 
ation of the slipping micro-j unction during the time it 
takes for the interface to cool down [43]. The effects of 
such temperature rises on friction have recently received 
renewed attention (see e.g. [62]), but remain poorly un¬ 
derstood. In an attempt to include such thermal pro¬ 
cesses in our model, we recognize that they will lead to 
time- rather than distance-controlled relaxations, which 
will distribute the shear force drop in time. In order for 
the interface to continue bearing the normal forces ap¬ 
plied to it, the micro-j unctions that relax are replaced by 
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new, pinned junctions that bear a small tangential force 

/new 


B. Technical aspects 

The physical aspects described above have been mod¬ 
eled in a simple way using the following assumptions. We 
consider the rough frictional interface between a horizon¬ 
tal track and a thin linear elastic slider (Fig. la). The 
slider has mass M and sizes L and H in the horizontal 
{x) and vertical (z) directions, respectively. We present 
the values of all parameters in Table I in Appendix C. 
The bulk elastodynamics of the slider are solved using 
a square lattice of blocks connected by internal springs 
(Fig. lb) [9, 58]. The slider is divided into a square lat¬ 
tice of N = NxNz blocks of mass m = M/N. Blocks are 
coupled to their four nearest neighbours and their four 
next-nearest neighbours by springs of equilibrium lengths 
I = L/{Nx — 1 ) = H/{Nz — 1 ) and a/ 2 / and stiffnesses 
k and /c/2, respectively, giving an isotropic elastic model 
with Poisson’s ratio 1/3. The force exerted on block i 
by block j is thus kij{rij — Uj)^^ when blocks are con¬ 
nected, 0 otherwise, where x = = Xj — x^, 

Vij = |Ax^j| and kij and lij are the stiffness and equilib¬ 
rium length of the spring connecting blocks i and j. Block 
oscillations are damped by introducing a viscous force 
r]{-kj—ici) on the relative motion of connected blocks. We 
chose the coefficient 77 = \/0.1km so that blocks are un¬ 
derdamped and event-triggered oscillations die out well 
before the next event. The non-frictional boundary con¬ 
ditions are the same as setup 2 in [9], which was also 
used in [27]. The top blocks are submitted to uniformly 
distributed, time-independent vertical forces The 

bottom blocks lie on an elastic foundation of modulus 
kf = /c/2, i.e. each block is submitted to a vertical force 
of amplitude pi = kf \zi\ if Zi < 0 or 0 otherwise, where 
Zi is the vertical displacement of block i. Both verti¬ 
cal boundaries are free, except for a horizontal driving 
force Ft = K{Vt — x^) applied on the left-side block 
situated at height h above the interface, where Xh is the 
^-displacement of this block. This models a pushing de¬ 
vice of stiffness K driven at a small constant velocity V. 

The multi-contact nature of the interface is modeled 
through an array of Ng tangential springs represent¬ 
ing individual micro-junctions, attached in parallel to 
each interfacial block (Fig. Ic) [25, 26]. The individ¬ 
ual spring behavior is as follows (Fig. If, [27, 47]). A 
spring pinned to the track stretches linearly elastically 
as the block moves, acting with a tangential force /t 
on the block. When the force reaches the static fric¬ 
tion threshold /thres (we neglect aging, so that /thres is 
time independent), the micro-junction ruptures and the 
spring becomes a slipping spring acting with a dynamic 
friction force /t = /slip- After a random time Ir drawn 
from a distribution T{tR), the slipping spring relaxes. 
It is replaced immediately by a pinned, unloaded spring 


(/new = 0) and a new cycle starts. Here we use T{tR) as a 
simplified way of modeling the distribution of times after 
which micro-j unctions relax. Due to the variety and the 
complexity of the underlying thermal processes, we did 
not try to derive T{tR) for a specific situation. Rather, 
we chose to model T{tR) in the simplest way, as a Gaus¬ 
sian with average time (tR) and width SIr. The shape 
of T(tR) is not crucial: we obtain qualitatively similar 
results with an exponential distribution. The width of 
T{tR) is the only source of randomness in our model and 
causes the interface springs of a block to evolve differently 
from each other. 

The 2N equations of motion are solved simultaneously 
using a leapfrog / velocity Verlet integrator [76] on a 
uniform temporal grid of resolution At. 


C. Relationship to other models 

1. Bulk modeling 

A spring-block discretisation of the bulk elasticity is 
particularly convenient for models where the friction is 
described as an ensemble of micro-junctions rather than a 
continuum law, because the blocks provide natural units 
on which to couple the frictional and the bulk elastic 
behaviour. 

Let us note that like finite element (FEM) and finite 
difference (FDM) methods, the spring-block discretiza¬ 
tion satisfies the equations of linear elasticity. In par¬ 
ticular, longitudinal (P) and shear (S) waves in the bulk 
propagate with the correct speeds and the right reflection 
and refraction properties [58]. To verify our implemen¬ 
tation we checked that the code reproduces the expected 
bulk wave speeds. 

The choice of the spatial resolution, the size of a single 
block, is made according to the following physical argu¬ 
ments. First, as discussed for example by Persson [44], by 
Carol! and Nozieres [63] and by Braun et al. [64], below a 
characteristic length scale A, called the elastic screening 
length, the interface behaves rigidly. A is thus the max¬ 
imum block size allowing for a correct representation of 
the elasticity of the interface. For a linear elastic rough 
interface, A ^ d^/a, with a the typical lateral size of mi¬ 
crocontacts and d the typical distance between them. For 
micrometer-ranged roughnesses, we expect a ^ 1 pm and 
d ^ 10 — 100 pm, yielding A ^ 0.1 — 10 mm. Second, 
the frictional behaviour of each block is then a statisti¬ 
cal average over the many micro-junctions connected to 
it. This statistical approach is increasingly relevant for 
larger blocks involving more junctions. One thus looks 
for the largest possible block size. Combining both re¬ 
quirements, A appears as the natural block size for such 
spring-block models. 
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2. Interface modeling 


The rate-and-state picture of friction, which includes 
displacement-controlled disappearance of micro-contacts, 
has proved to be adequate for slow (typically up to 
100 pm/s range) sliding in a variety of materials. Such 
slow velocities imply a negligible temperature rise of the 
interface and thus a slipping state the duration of which 
is of purely geometrical origin, i.e. it is controlled by a 
length scale of the order of the mean micro-contact size. 

Here, we focus on a drastically different situation, in 
which the transition from static to kinetic friction is ex¬ 
tremely short (millisecond range) and is accompanied by 
fast slip (100 mm/s range). As recognized in recent ex¬ 
periments by Ben-David et al. [59, 65], the sudden rup¬ 
ture of the interface and its subsequent slip will generate 
a significant heating of the interface, sufficient to melt 
the broken micro-asperities. In these severe conditions, 
precise knowledge about the micro-contact behavior at 
the millisecond time scale is currently lacking. 




t/p 

— t^ 


Time (s) 


0 



In our description, we acknowledge the fact that the 
onset of sliding is far from the slow steady sliding sit¬ 
uation. During the time in which the interface is sig¬ 
nificantly heated, thermal activations for the transitions 
from slipping to pinned states are highly probable. Such 
activation is classically described by time-rates, rather 
than by displacement-related quantities. It is therefore 
natural to propose an alternative picture that incorpo¬ 
rates the possibility that transitions between the slipping 
and the pinned states can be controlled by a time. 


FIG. 1: Color) Sketch and behaviour of the multiscale 
model, (a) Slider and external loading conditions, (b) 
Spring-block network modeling elastodynamics. (c) 
Surface springs modeling friction on a block, (d) 
Macroscopic loading curve, the ratio Ft/Fn of driving 
shear force to total normal force, (e) Mesoscopic loading 
curve, the ratio r/p of shear to normal stress on a block, 
(f) Microscopic friction model for the spring loading 
curve, the ratio fr/I n of friction to normal force for 
one spring (/tv =p/Ns). Figure adapted from [27]. 


We emphasize that experimental data did show that 
the dynamics at the onset of sliding involves a transi¬ 
tion from fast to slow slip which occurs after a constant 
time, rather than a constant displacement [see 59]. This 
timescale controls the fast dynamics with which the in¬ 
terface comes back to a fully pinned state after slipping. 
It thus drastically differs from the classical time scale for 
aging, also found in [59], and which controls the slow 
strength recovery of the interface when it is at rest. 

Note, however, that we would not be surprised if a 
length scale would also be relevant to the dynamics of 
the onset of sliding in the same experimental system. 
However, the study of a complete model involving both 
a time scale and a length scale is far beyond the scope of 
the present work. 

In terms of modeling approach, let us also stress the 
fact that a number of reference models from the litera¬ 
ture considered, before us, time-controlled transitions be¬ 
tween micro-junction states [see e.g. 25, 26, 43, 47, 66]. 
As shown in Thpgersen et al. [47], our friction law is 
actually one particular case of a more general family of 
models. 


III. LOADING CURVES AND SLIP DYNAMICS 
A. Loading curves 

The evolution in time of the driving force Ft is called 
the loading curve. It is readily measured in experiments, 
and is used to characterize the motion as smooth slid¬ 
ing, or regular or chaotic stick-slip. The loading curve 
for our chosen set of parameters is shown in Fig. Id. It 
starts with an initial buildup from zero load. This ini¬ 
tial buildup is linear, as we apply the load through a 
linearly elastic spring, one end of which is pushing the 
sample while its other end is being driven with a con¬ 
stant speed. From about 0.2 s the linear increase is in¬ 
terrupted by small drops in the driving force. These are 
associated with precursors: rupture events confined to 
only part of the interface [see e.g. 6-10, 12, 14, 25, 67]. 
The first event where the entire interface breaks and the 
slider moves macroscopically occurs at about 0.45 s and 
is seen in the loading curve as a larger drop in Ft whose 
duration is not resolved on this figure. Then the system 
enters regular stick-slip, with alternating events showing 
partial and full breaking of the interface (small and large 
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drops in Ft). In the model, as in experiments, the be¬ 
haviour on larger scales emerges from the interactions on 
the scales below. Figure If illustrates the single junction 
law described in Section 11. 

Each mesoblock is coupled to the track with many 
junctions (here we use Ns = 100). An example of the evo¬ 
lution of the net force in the junctions, which apart from 
a small inertial term is equivalent to the sum of forces on 
the block from its neighbours, is shown in Fig. le. The 
data is from the stick-slip phase. The temporal structure 
is rich in detail, but can overall be described as periods 
of very slow increase (the driving is slow compared to 
the internal dynamics) separated by faster increase and 
decrease when rupture fronts pass. 

The force maxima attained on the mesoblock level are 
smaller than the sum of the individual junction thresh¬ 
olds, which would give rjp = 0.4. The explanation is 
simple: due to the disorder in the individual stretching 
states, some junctions will break before others, and so 
they will not all contribute their maximum force simul¬ 
taneously (see detailed discussion in [27, 47]). Similarly, 
the maxima in the macroscopic loading curve of Fig. Id 
are smaller than those on the mesoblock level, as the 
mesoblocks do not reach their individual maxima at the 
same time. This difference between macroscopic and lo¬ 
cal friction coefficients has been discussed in e.g. [7- 
9, 26, 33, 68]. That the overall strength in a system is 
usually smaller than the sum of the individual strengths 
of the constituents is familiar from other fields, for ex¬ 
ample fracture mechanics [see e.g. 69] and fiber bundle 
theory [see e.g. 70]. 


B. Block slip dynamics and a slow slip mechanism 

Figure 2a shows the slip dynamics of a block in a par¬ 
tial slip event and a block in a full sliding event. We 
find both fast and slow slip regimes of the motion, in ex¬ 
cellent agreement with the experiments reported in [59]. 
The initial fast slip regime begins after the passage of a 
fast rupture front. In both these cases the fast slip is 
followed by the block coming nearly to rest (no visible 
increase in net slip between 450 and 550 ps) and then by 
a slow slip regime with roughly linear increase of slip vs 
time, i.e. constant slip speed. The slow slip regime can 
end in two ways. In full sliding events, slow slip changes 
back to fast slip when the slider enters the full sliding 
regime. For arresting events, the slow slip regime ends 
when the block comes to rest. 

The initial fast slip regime corresponds to an inertial 
motion of the block when a large number of junctions 
break in a short time interval as the rupture front passes 
by. The net friction force is rapidly reduced, bringing 
the block out of mechanical equilibrium. The result is 
a large positive acceleration (in the direction of the net 
force due to the neighboring blocks). The inertial nature 
of this motion is demonstrated in Fig. 6. 

The subsequent slow slip observed in the model has a 




FIG. 2: Color) (a) Example slip profiles from a partial 
slip and a full sliding event. Line without markers 
(black): slip profile for block at x = 0.161/ with 
to = 0.7872 s, i.e. the partial slip event preceding event 
III. Line with markers (grey/multicolored): slip profile 
for block at X = 0.341/ with to = 1.0571 s, i.e. event II 
(more details of this event in Fig. 7). In both cases we 
chose a block located near the middle of the region 
where the rupture front speed was fast. For blocks 
closer to the fast-slow transition or to the front arrest 
point, the amplitude of fast slip is smaller, (b) Sketch 
showing how the relaxation of force associated with the 
junctions’ relaxation from the slipping (s) to the pinned 
state (p) can lead to a slow slip motion of the block. 


different physical origin that was explained in [27, 47]. 
This slow slip mechanism is illustrated in Fig. 2b. When 
/new < /slip the uct frictiou force on a block is reduced 
slightly whenever a spring leaves the slipping state, yield¬ 
ing a small positive acceleration as the net friction drops 
below the net external force from neighboring blocks. 
The friction reduction is soon balanced by the changes in 
the pinned junctions and the external forces on the block 
as it slowly moves. This slow slip mechanism is present 
as long as some junctions are going from the slipping to 
the pinned state and /new < /slip, but it is masked by the 
fast slip while the fast slip lasts. The dependence of slow 
slip speed on model parameters was discussed in detail 
in [27]. 


IV. FRONT TYPE RESULTS 

In this section we present our results on the qualitative 
features of rupture fronts. That is, we discuss the condi¬ 
tions under which we observe fast rupture, slow rupture, 
and the transitions between these regimes. 
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ftaction of pinned springs 


100 % 



(c) 



FIG. 3: Color) Two interface-sized events, (a) A 
fast-slow-fast event (I in Fig. Id). Spatiotemporal plot 
of the fraction of pinned springs, (b) A fast-only event 
(III in Fig. Id) shown as in (a), (c) Rupture front speed 
Vc vs. front location for both events. Block rupture is 
defined to occur when 70% of interface junctions have 
broken (white dashed line in the colorbar). Front speed 
is measured as the inverse slope of the rupture line 
(indicated by arrows in (a) and (b)) using the endpoints 
in a five-point-wide moving stencil. Figure adapted 
from [27]. 


A. Rupture front characteristics 

With the driving force applied on the trailing edge of 
the slider, the blocks near the trailing edge are the first 
to reach their effective static friction thresholds. A block 
that slips increases the load on its neighbors, which can 
start to slip in turn. The rupture front tip, i.e. the 
boundary between a region of stuck blocks and a region 
of slipping blocks, then propagates away from the nucle- 
ation point. If the rupture arrests before reaching the 
leading edge of the slider the event is called a partial slip 
event; if the leading edge is reached, we use the name 
full/global sliding event; precursors are partial slip events 
that occur before the first global sliding event. 

In the interaction law the distinction between pinned 
and slipping states is made on the junction rather than 
the block level. We therefore define slipping on the block 


level to mean that a certain fraction of the block’s junc¬ 
tions are in the slipping state. This criterion is robust 
to the choice of threshold fraction, because as is seen in 
Fig. 3 and other figures (7ab, 8, 21), it is typical for a 
block to go from having nearly all its junctions pinned 
to having nearly all of them broken in a time short com¬ 
pared to the other time scales in the simulation. For 
the events in Fig. 3 the time to go from 80% of springs 
pinned to 20% of springs pinned is approximately 0.03 ms 
in the fast part of the fronts and 0.3 ms in the slow part 
of the fronts. We have used a threshold value of 30%, 
that is, the start time of block slipping is taken as the 
instant when the fraction of pinned junctions dropped 
below 30%. 

Figure 3a shows a rupture front whose speed Vc changes 
from fast {vc ^ Cs/3) to slow {vc ^ Cg/lOO) and back to 
fast again. The left-travelling front that starts when the 
primary front is reflected from the leading edge re-breaks 
the junctions that had healed behind the front tip. Fig¬ 
ure 3b shows a rupture front that is fast across the entire 
interface. By defining the transition to block sliding as 
above, the location of the front tip in time can be mea¬ 
sured. The local front speed is then the ratio of the 
distance travelled by the front to the time for that prop¬ 
agation. Because of disorder in the junction state along 
the interface remaining from earlier events, the propa¬ 
gation time from one block to the next can vary signifi¬ 
cantly. We have found that using only the end-points in 
a 5 blocks wide moving stencil (Appendix A) gives the 
best balance between robustness and spatial resolution 
in the calculated front speed. The results are shown in 
Fig. 3c. 


B. The influence of front type on the loading curve 

The two events in Fig. 3 are different. One is a fast- 
slow-fast front, the other is fast across the entire inter¬ 
face. Both events are marked in Fig. Id, and the associ¬ 
ated drops in Ft are seen to have approximately the same 
amplitude. This indicates that the details of the front 
propagation do not influence the loading curve strongly, 
at least not for these events. Nevertheless, the slow front 
propagation does have a signature in the loading curve 
that appears when we zoom in on a few events as in 
Fig. 4. Namely, the drop in the loading force has a sig¬ 
nificant change in slope while the slow front lasts, which 
distinguishes it visually from the force drop associated 
with a fast-only event. 

Let us consider the evolution of Ft in more detail. 
As long as the slider remains pinned. Ft increases with 
the motion of the driving stage (the driving stage moves 
the end of the driving spring that is not attached to 
the slider). Ft decreases only when the point on the 
slider where the driving spring attaches, the trailing edge, 
moves away from the driving stage. This occurs in two 
distinct ways. First, the trailing edge moves away from 
the driving stage when the slider deforms in compression 











during the passage of a rupture front, which happens 
both for partial slip and full sliding events. Second, the 
trailing edge moves with the rest of the slider when the 
entire interface is slipping in a full sliding event. For full 
sliding events we define the boundary between rupture 
front passage and full sliding as the moment the rupture 
front reaches the leading edge, seen e.g. in Fig. 3. The 
relative amplitude of the Ft reductions associated with 
each of these two motions of the trailing edge depends 
on the relative stiffnesses of the slider and the driving 
spring. When the slider is stiff compared to the driving 
spring, so that little motion of the trailing edge can oc¬ 
cur unless the entire slider moves, the drop in Ft is only 
appreciable during the sliding part of full sliding events. 
With a softer slider, the deformation occurring in par¬ 
tial slip events and during the rupture front passage in 
a full sliding event accounts for a larger fraction of the 
net reduction in Ft. As seen by the relative amplitudes 
of the force drops in Fig. 4, with the present parameters 
the (trailing edge) slip associated with the slider defor¬ 
mation accounts for about a fifth of the net slip of full 
sliding events. 

Further quantification of this feature is presented in 
Fig. 5, where the force drops associated with either the 
entire event or with the rupture front passage only are 
presented for all events during the developed stick-slip 
regime. We find that the loading force drop occurring 
during the rupture front passage, regardless of the type 
of event, is always much smaller than the one associated 
with sliding. The net force drop in full sliding events 
has comparable amplitudes for fast-slow-fast and fast- 
only events. This is because the block motion during 
full sliding accounts for the largest part of the net block 
motion in an event, and after the rupture front reaches 
the leading edge, the distinction between fast-slow-fast 
events and fast-only events is unimportant. 


C. Fast slip and fast front speeds are inertial 

In Fig. 2a we showed example block slip profiles and 
the regimes of fast slip, slow slip and full sliding. In sec¬ 
tion IIIB we focused on the slow slip part, and discussed 
how the junction evolution law leads to a slow slip mech¬ 
anism on the block level. In Fig. 3 we showed that the 
model exhibits both fast-only and fast-slow-fast fronts, 
and explained how we measure the front speed. In this 
section we demonstrate that both the fast slip part of 
the block slip evolution and the fast front speed are of 
inertial origin. That is, like the bulk wave speeds, these 
speeds scale as where p = Mj is the mass 

density of the system. 

To isolate the effect of inertia from the stress and fric¬ 
tional state at the interface we have performed four sim¬ 
ulations starting from the same state (Appendix C), but 
with p decreased to 1, 1/2, 1/4 and 1/8 of its value in 
Table I; we changed the mass and kept the system size 
constant. Figure 6a shows the fast slip dynamics for four 



FIG. 4: A zoom on the macroscopic loading curve 
Ft/F^. The drop in a full sliding event is comprised of 
two parts: first, the drop associated with deformation of 
the slider during rupture front passage (dashed arrow); 
second, the drop associated with motion of the entire 
slider (full arrow minus dashed arrow). The amplitude 
of the first drop is approximately the same as the 
amplitude of the drop associated with a partial slip 
event (see Fig. 5). 


neighboring blocks located within the part of the inter¬ 
face where the front speed was high as the front passed, 
for each of the four simulations. The figure clearly shows 
that the fast slip speed was modified by the change of 
density. Figure 6b demonstrates that a scaling col¬ 

lapses the data. The reference time trup for each block 
is the time of block rupture as defined in Fig. 3, and the 
block slip is measured with respect to the block position 

Gup- 

Figure 6c shows the front speed as a function of po¬ 
sition along the interface for the four simulations, while 
Fig. 6d shows the same data with the front speed rescaled 
by Xhe fast front speeds are collapsed onto each 

other by this rescaling. The slow front speeds, in con¬ 
trast, remain the same in all four simulations (Fig. 6c), 
and are split from each other by the scaling. This indi¬ 
cates a non-inertial origin of the slow front speed, which 
is the topic of the next section. 


D. Prom slow slip to slow fronts 

In the model the propagation of fast rupture and slow 
rupture is governed by different mechanisms. To demon¬ 
strate this we start from the fast-slow-fast event in 
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FIG. 5: Color) The force drops in the loading curve 
(Fig. Id) grouped according to event type. We include 
partial slip events and full sliding events occurring 
between t = 0.63 s and t = 1.50 s (the developed 
stick-slip regime). The drops occurring during the 
passage of the rupture front have comparable amplitude 
between the partial slip events, fast-slow-fast events 
and fast-only events. Also, the net force drop in full 
sliding events have comparable amplitudes for 
fast-slow-fast and fast-only events. 


Fig. 7a. We then modify the junction law by increas¬ 
ing /new to /slip, so that there is no change in the friction 
force when a slipping junction relaxes and is replaced by 
a pinned junction. This does not affect the initial force 
relaxation of pinned junctions reaching /thres? but it does 
affect the second relaxation of slipping junctions relaxing 
back to the pinned state. In practice, this change turns 
off the slow slip mechanism that was discussed in Sec¬ 
tion IIIB. The result of restarting from the same state as 
in Fig. 7a and with only this modification to the junction 
law is seen in Fig. 7b. The first, fast part of the rupture 
is unaffected, but the slow rupture is suppressed, and the 
front stops where the fast-slow transition used to occur. 

Figure 7c illustrates how slow slip in the region behind 
the front tip can result in additional front propagation: 
as the blocks at and behind the front tip move towards 
the stuck region, the external forces on the block just 
ahead of the front tip increase. Two outcomes are possi¬ 
ble. Either the block remains stuck and the front arrests, 
or the forces on this block eventually overcome its effec¬ 
tive static friction threshold, and the block starts to slip, 
moving the front tip one block ahead. This last part of 
the explanation would be the same for fast slip and fast 
fronts. The difference between fast and slow fronts is 


a 


m 


x10® -- x10' 



Time, t — (s) (t — trup)/\/p (arb. units) 




FIG. 6: Golor) The fast slip and the fast front speeds 
scale with inertia, (a) Block slip motion for four 
neighboring blocks within the fast front region, for four 
simulations of the same fast-slow-fast event. In each 
simulation the initial state is the same, but the mass 
density p is different between the simulations. The 
block slip is measured from the rupture time trup of 
each block as defined in Fig. 3. Lines with the same 
color and marker are from the same simulation. Lines 
with the same line style represent the same block in 
different simulations, (b) Rescaling the time of slip with 
collapses the data in (a), (c) The rupture front 
speed as a function of position for the same simulations 
as the data in (a) (corresponding colors and markers). 
The change of density modified the fast front speed, 
while the slow front speed remained nearly unchanged. 

(d) Rescaling the front speed data in (c) by 
collapses the fast front speed measurements, but splits 
the slow front measurements. 


thus simply a matter of the origin of the underlying slip 
motion, which is inertial for fast fronts (see Fig. 6), and 
is related to the intrinsic arrest dynamics of the interface 
for slow fronts (see [27, 47]). 

Our understanding of the transition from fast to slow 
fronts is that the fast propagation stops at the same point 
in both Fig. 7a and b. In Fig. 7a, slow slip becomes 
important when the fast slip ends, and results in a slow 
front propagating. In Fig. 7b, the slow slip mechanism is 
turned off, and the event is over once the initial fast slip 
ends. 

We emphasize that the physical origin of the slip on 
the block level is unimportant to the way slip behind the 
front tip leads to rupture front propagation. We showed 
in [27] that a mechanism for slow slip completely different 
from the one in the present model would also lead to slow 
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FIG. 7: Color) (a) and (b) Event II restarted at 1.0571 s 
with driving speed 17 = 0, shown as in Fig. 3. Slow slip 
is either (a) active or (b) turned off. (c) A sketch 
showing how slip behind the front tip leads to front tip 
propagation. Panels (a) and (b) adapted from [27]. 


mind, the strategy for turning fast-slow-fast events into 
fast-only fronts and the other way around was to mod¬ 
ify the width of the initial force distribution of selected 
blocks along the interface. 

Slow slip becomes important to front propagation only 
after a fast front stops. We can keep the initial fast part 
of the event in the re-simulation equal to the original 
event by leaving the stress state unchanged. By increas¬ 
ing the width of the junction force distribution we make 
the interface weaker, which may enable fast propagation 
across the whole interface. An example of this is shown 
in Fig. 8b. 

Conversely, in Fig. 8d, by significantly decreasing the 
width of the initial junction force distribution, we made 
the interface stronger along a part of the interface, so that 
the original fast-only propagation was stopped. This al¬ 
lowed the slow slip mechanism to become important, and 
the front was turned into a fast-slow-fast one. 

We note that in both cases, the stress state in the 
original and modified simulations are the same, and the 
only change is in the width of 0 (/t) in the region where 
it is modified. To keep the initiation of the events the 
same in the modified and original simulations, we did 
not modify (j) near the trailing/left edge of the system. 


front propagation. 


E. Junction force distribution affects front type 
selection 

We argued in the previous section that even when the 
slow slip mechanism is active, its effect on the rupture 
front speed can be masked by fast slip and the associated 
fast front propagation. To illustrate this (Fig. 8), we will 
here re-simulate a fast-slow-fast event with the initial 
state modified so that the event becomes fast-only. Simi¬ 
larly, we will re-simulate a fast-only event with the initial 
state modified so that the event becomes fast-slow-fast. 

To understand these results, recall the connection be¬ 
tween the distribution 0 (/t) of forces among the junc¬ 
tions attached to one block and the corresponding ef¬ 
fective static friction coefficient of this same block. As 
shown numerically in e.g. [27] and theoretically in e.g. 
[42, 47] and Appendix B, for friction, and more generally 
for the rupture of heterogeneous systems in which a num¬ 
ber of junctions are loaded in parallel [70], the maximum 
load that an interface can bear is related to the width 
of the load distribution and/or threshold distribution of 
the various junctions. Homogeneous systems (vanishing 
distribution width) have the maximum possible macro¬ 
scopic rupture threshold because all junctions will con¬ 
tribute with their maximum force when collective rupture 
is reached. In contrast, in heterogeneous systems (finite 
width) the weaker and/or initially more highly loaded 
junctions will break first, so that when macroscopic rup¬ 
ture occurs, only a fraction of the initial population of 
junctions will contribute to the total force. With this in 


F. Front type phase diagram and its predictive 
power 

We have seen that increasing (decreasing) the width of 
the junction force distribution makes the interface weaker 
(stronger) and that this favours fast (slow) front propaga¬ 
tion. It also makes intuitive sense that higher prestress, 
hence smaller distance to the breaking threshold, would 
favor fast front propagation (we return to this in more de¬ 
tail in Section V). We have performed simulations where 
these two parameters are varied systematically (see Sec¬ 
tion VB and Appendix C), and the observed front types 
are presented in Fig. 9. In the arresting region, the fronts 
are partial slip events, that is, they stop before reaching 
the leading edge (some of them stop early, some almost 
reach the leading edge). The region labelled slow includes 
all those events that have a slow front part, even if the 
event is fast along most of the interface. They share the 
characteristic that the events would arrest in the absence 
of the slow slip mechanism. In the fast region, events are 
fast across the entire interface. 

In Section V we discuss the transient behavior of the 
rupture fronts. For now, let us stress that while the ar¬ 
resting, slow and fast regions of the rupture fronts are 
robust in their relative positions (the fast front region 
is found at higher values of normalized prestress tq and 
junction distribution width a than the slow front region, 
which is itself found at higher values than the arrest re¬ 
gion), the precise locations of the boundaries between 
them depend also on how the events are triggered. For 
example, it is possible to prepare two simulations with 
the same stresses and distributions in the propagation 
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FIG. 8: Color) Junction force distributions affect front 
propagation. Main panels: events shown as in Fig. 3. 
Insets: histograms of junction force distributions, (a, b): 
A fast-slow-fast event turned into a fast-only event by 
increasing the width of the stretching distribution while 
maintaining the shear stress (the modified distributions 
are uniform). The insets are taken at t = 1.054 s. (c, d): 
A fast-only event turned into a fast-slow-fast event by 
decreasing the width of the stretching distribution along 
part of the interface. The shear stress was unchanged 
and the new distribution was uniform with the width 
shown in the inset of panel (d). The insets are taken at 
t = 0.824 s. All insets: for each block along the 
interface, a color coded histogram of 0 (/t)- The 
vertical axis shows the force level in individual springs, 
which extends up to /thres- Upper white line: the level 
/thres, which is different for each block because it varies 
with normal force. Lower white line: /t = 0 N. The 
color scale denotes the fraction of springs found at each 
value of /t using an arbitrary bin width. Offset data: 
fraction of slipping springs at 1.054 s (a, b), at 0.824 s 
(c, d). Panels (a) and (b) adapted from [27]. 



FIG. 9: Golor) A front type ’’phase diagram” based on 
simulation of homogeneous prepared interfaces. Front 
type observed vs initial width a of junction force 
distribution and prestress 

To = (ro/p- /siip//jv)/(/thres//iv “ /siip//iv)- Initial 
states in the forbidden region would have had some 
junctions stretched beyond their breaking threshold and 
are therefore excluded. Figure adapted from [27]. 


region, but with different stresses in the triggering region 
(these regions are defined in Fig. 12c and in Appendix G) 
such that the simulation with the more highly stressed 
triggering region produces a fast-only front and the other 
a fast-slow-fast front. For this reason, all the simulations 
in Fig. 9 have the same settings in the triggering region. 

We also wish to note that the front type phase diagram 
is not a local measure. For example, as mentioned above, 
many of the events that populate the slow front region are 
of the fast-slow-fast type, meaning that even though the 
prestress and junction force distributions are the same 
for all the blocks in the propagation region, the rupture 
passes some of them as a slow front and others as a fast 
front. However, if these limitations are kept in mind, 
the diagram is still a powerful tool for guiding our intu¬ 
ition. For example, in an interface where there is a sud¬ 
den change of initial conditions within the propagation 
region, the phase diagram tells us what change to expect 
in the rupture front. In Fig. 10 we show three simula¬ 
tions. The reference simulation (Fig. 10a) for this figure 
is the simulation behind the point at (0.088, 0.0875) in 
Fig. 9. The two other simulations show that when there 
is a sudden increase of prestress or junction force distri¬ 
bution width in the propagation region, this can lead to 
a transition from slow to fast front propagation right at 
the location where the change occurs. Gonversely, some 
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FIG. 10: Color) Simulations showing that even though 
the phase diagram in Fig. 9 is an inter face-wide and not 
a local measure, the change in front characteristics 
upon local changes to prestress and spring stretching 
width are consistent with the predictions of the phase 
diagram, (a) The original simulation with homogeneous 
stress and junction distribution states in the 
propagation region, (b) A step increase in the initial 
shear stress within the propagation region leads to a 
slow-fast transition, (c) A step increase in the width of 
the junction force distribution within the propagation 
region also leads to a slow-fast transition. Top: events 
shown as in Fig. 3. Middle: Spatial distributions of 
initial prestress tq/p (before event triggering) is shown 
with drawn line. The grey region is the triggering 
region (see Appendix C). Bottom: histograms of 
junction force distributions shown as in Fig. 8, taken 
after the triggering of the event, at t = 0.0207 s. 


(but not all) simulations with a change of opposite sign 
(from a high to a low prestress, for example), lead to a 
fast-slow transition. 


V. FRONT SPEED RESULTS 

The rupture front speed is the speed at which the slip¬ 
ping region grows into (invades) the region which is still 
stuck. In the previous section we distinguished slow rup¬ 
ture and fast rupture and found that they are governed 
by different mechanisms. In this section we go on to con¬ 
sider the difference in speed between two events that are 
either both fast, or both slow. 

In the language of fracture mechanics, rupture occurs 
when the energy available at the front tip reaches the 
energy required to break the contacts there. Thus, the 


speed at which rupture propagates depends on the en¬ 
ergy that is already present (related to the stress state), 
the energy level needed to break the contacts (related to 
the local strength) and on how quickly the missing en¬ 
ergy can be supplied (which depends on the slip motion 
behind the tip and is therefore transient). Although it 
has been shown that fracture mechanics satisfactorily de¬ 
scribes fast shear rupture events [14, 22], in the model, it 
is more straightforward to argue in terms of forces rather 
than energies; however, we will see that considerations 
similar to those in fracture mechanics apply. 


A. Front speed is transient 

If the front speed did not have transients, or if the tran¬ 
sients were short compared to the length scale at which 
we study the front propagation, the rupture front speed 
would be uniquely determined by the local state of the 
interface at the rupture tip, i.e. local stresses, strength, 
stiffnesses etc. Indeed, various local interfacial param¬ 
eters have been shown to be correlated with the local 
front speed (see e.g. [9, 16, 32]). However, our simula¬ 
tions show that the speed at any point depends not only 
on the state at that point, but also on the region the 
front has just passed through, that is, the front speed 
has transients. Near the trailing edge, within the length 
required for the speed of the newly nucleated front to 
converge, the region behind the front that influences the 
propagation extends back to the region where the front 
was triggered, and so in addition to the state of the inter¬ 
face, the precise way rupture starts also is part of what 
determines front speed. 

To illustrate the transient nature of the front speed. 
Fig. I la shows the front speed as a function of posi¬ 
tion for two simulations in which the local stress state 
and local frictional strength are homogeneous along the 
propagation region (see details about how these single 
event simulations were set up in Appendix C). We first 
consider the shorter of the two systems (magenta dotted 
line). We observe that the front speed changes through¬ 
out, that is, even though the local state is the same along 
the slider, the front speed is not. The length of this sys¬ 
tem is Nx = 57, Lx = 140 mm, and the front is in the 
transient part of propagation throughout. Note that this 
is a fast-only front; even though the front speed starts 
low, its propagation does not depend on the slow front 
mechanism. 

To investigate the convergence length of the transients 
and the spatial extent over which edge effects dominate 
we performed another simulation. This simulation differs 
only in the total length of the sample, which we increased 
by a factor of five. This gives the black line with circular 
markers in Fig. I la. We observe that the two rupture 
events have very similar front speeds. When the front 
tip approaches the leading edge, there is a change in cur¬ 
vature that we interpret as an edge effect. This interpre¬ 
tation is supported by the fact that the effect is present 
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near the leading edge of the sample for both the short 
and the long system, but the transient effect observed 
at the end of the short sample, around blocks 50-57, is 
not observed at blocks 50-57 in the longer sample, where 
these blocks are far from the edge. 

Figure 11b, where front speed is plotted against inverse 
position along the interface, shows that the front speed 
does converge to a finite value. Note that we do not ex¬ 
pect a 1/x behaviour to hold in general; rather, we plot 
front speed versus inverse position because it is the sim¬ 
plest function that brings very large x values close to the 
ordinate axis, and thus facilitates a crude extrapolation 
to the expected front speed reached in very long samples. 

Figure 11c and d also indicate convergence. They 
shows two events. The local stresses and strengths are 
the same in both, but the initialization of the events is 
different. This leads to front speeds that are initially also 
different, but which converge to the same value. 


B. Front speed depends on local stress state 

While we and others [32] argue that the ratio tq/p of 
shear to normal stress at the interface as the rupture front 
begins provides insufficient information for predicting the 
front speed on its own, it remains true within the model, 
as is also seen in experiments [16], that the front speed 
has a strong dependence on To/p. Higher prestress tq 
results in higher front propagation speed. This result 
agrees with our intuition: for a given strength, a region 
that is more highly prestressed requires less stress change 
to start slipping, and it releases more energy when doing 
so. Both favor higher front speed. This is also consistent 
with Fig. 9, which showed that increasing the initial shear 
stress in the propagation region will bring a system into 
a regime where fronts are fast across the entire interface. 

Figure 12a shows the front speed for events in which 
the initial stress state was varied systematically. A simi¬ 
lar variation in the junction force distributions’ width a 
is the topic of Section V C. To isolate the effect of varying 
prestress from other variations that arise in the full simu¬ 
lations we used the protocol illustrated in Fig. 12c, which 
shows the initial shear and normal stress states. The nor¬ 
mal stress p is homogeneous, i.e. taking the same value 
everywhere along the interface. The shear stress is also 
homogeneous within each of two regions. The region on 
the left is called the triggering region, because we start 
each event by simultaneously and artificially detaching 
all junctions for the blocks within this region. This en¬ 
sures that the starting conditions are the same for all the 
fronts. We always use the same shear stress level in this 
region, Tq = ( tq/p - /slip//Ar)/(/thres//Ar - /slip//iv) = 
0.3. In the rest of the interface, which we call the prop¬ 
agation region, the initial shear stress is also uniform, 
but not necessarily the same as in the triggering region. 
For clarity only two prestress values are shown in the 
propagation region of Fig. 12c. 

In agreement with our argument above we observe in 








2000 

1500 

1000 

500 

0 

2000 

1500 

1000 

500 

0 




FIG. 11: Color) Results for transient behavior of fast 
rupture fronts. Front speed vs position (a) and inverse 
position (b) for simulations where the local state is 
homogeneous along the propagation region. The 
magenta/grey line with dotted markers is for a system 
of our reference length, the black line with circular 
markers for a system five times longer, (c) and (d) 
Difference in front speed for two systems of the same 
length (five times reference length) and the same initial 
state. Black line with circular markers: the same data 
as the black lines in (a) and (b), which was triggered by 
simultaneously breaking all junctions for all blocks in 
the triggering region as explained in Appendix C. 
Green/grey line with dotted markers: front triggered by 
driving from the trailing edge, which modifies the stress 
in the loading region prior to rupture. In all panels, 
rupture velocity is defined as in Fig. 3. In (b) and (d), 
the speed for the ten blocks closest to the leading edge 
are excluded from the plot, to avoid the region where 
the edge effects dominate. 


Fig. 12a that higher prestress (dark red line) corresponds 
to higher mean front velocities. 

Note that in the model the longitudinal wave speed 
Cl = 1677 m/s does not set an upper limit for the front 
propagation speed. This is consistent with the results 
of other experimental (see e.g. [71]) and numerical (see 
e.g. [32]) studies of bi-material interfaces, in which the 
maximum front speed is expected to be the one of the 
stiffest material. Here the track is perfectly rigid in the 
front propagation direction, so that we do not expect any 
conceptual limitation to the front speed. 















14 




Position, X (m) Prestress, tq 


0.4 

^ 0.3 

I? 

^ 0.2 
lo. 


/thres/ fi 


N- 


/siip/ In 

Triggering Propagation 
region region 


(c) 


0 0.02 0.04 0.06 0.08 0.1 

Position along interface, x (m) 


0.12 0.14 


FIG. 12: Color) Rupture front speed depends on local 
stress state, (a) Rupture front speed vs position along 
the interface for a series of simulations with various 
homogeneous initial shear stresses in the propagation 
region. The junction distribution width was zero, so 
these simulations would plot on top of the vertical axis 
in Fig. 9. They are all of the fast-only type. The 
transient nature of the fronts discussed in Section V A 
can be seen here as well, in the increase in front speed 
towards the right of the sample, (b) The front velocity 
data in (a) vs the prestress value for each simulation. 

Dots: mean values. Bars: min/max values. Each 
dot+bar corresponds to the line in (a) with the same 
color / grey scale + vertical extent. The prestress 
values are 0.22 to 0.94 in steps of 0.06. (c) Spatial 
distributions of prestress ro/p. More details in the text 
and in Appendix C. For these simulations, a time step 
of At = 5 • 10“^ s was used to obtain smooth 
measurements at high front speeds. Panel (c) adapted 
from [27]. 


C. Front speed depends on local strength 

In the same way as for the influence of stress state 
on front speed discussed above, it makes intuitive sense 
that for a given prestress, the front propagates faster if 
the frictional strength, i.e. the threshold stress for slip 
inception, is lower. In the model, the frictional strength 
is controlled to first order by the values of /thres and 
/slip oa the junction level, but even when these are the 
same everywhere, the local strength can vary. We show 


in Fig. 20 (Appendix B) how increasing the width, cr, of 
the initial junction force distribution, 0 (/t), results in 
a lower effective static friction threshold, and we expect 
that if conditions are otherwise the same, this will result 
in faster front propagation. In the case (not studied here) 
of a more complete model, variations in local strength 
would also be due to variations in individual strength of 
the junctions attached to each block. 

We observe in Fig. I3a that for the same prestress, 
higher a results in higher front speed. Within our sam¬ 
ple size, which is the same as in previous sections and was 
selected from experimental parameters, the front speed 
is changing throughout the front propagation. In order 
to investigate the effect of a quantitatively and to isolate 
it from the effect of prestress, we compare each event 
to a reference event having the same prestress state and 
with width a = 0. We measure along the entire inter¬ 
face and find average and min/max values. The results 
are presented in Fig. I3b. By rescaling the front speeds 
with the reference speed for a = 0 the similarity between 
the two series for different prestress becomes apparent. 
The range of results for each value of a due to the tran¬ 
sient nature of the front speed is large, but there is a 
clear trend that higher a (weaker interfaces) correspond 
to higher front speeds. The effect of a on speed can be 
very significant, giving up to a 50% relative increase for 
the data shown in Fig. I3b. 

Figure 13a and b show data for fast-only fronts. The 
influence of a on front speed is qualitatively similar for 
slow fronts. This is shown in Fig. I3c and d. For the 
fronts, which have both fast and slow propagation parts, 
we used visual inspection of the underlying data in Ap¬ 
pendix D to determine which part of the fronts to include. 
Our selection can be read off of the horizontal extent of 
the lines in Fig. I3c. Again, the effect of a on the front 
speed is significant, with an increase of about one order of 
magnitude when cr/(/thres ~ /slip) is increased from 0.12 
to 0.47. 


D. Front speed is proportional to slip speed 

The motion of the rupture front tip and the motion 
of the material of the slider are interrelated: the slider 
cannot move while the interface is in the pinned state (the 
slip depends on the front), and the rupture propagates 
as deformation of the slider transfers stress and energy 
to the front tip (the front depends on the slip). In [27] 
we demonstrated that the slow front speed in the model 
is proportional to the slow slip speed and worked out the 
constant of proportionality. In this section we show that 
the fast front speed in the model depends on the fast slip 
speed through the exact same relationship. 

Figure I4a shows the fast front speed vs fast slip speed 
for the blocks between x = 3 cm and x = 13 cm in a se¬ 
ries of simulations with controlled initial states and vari¬ 
ation in mass density, normal force, junction distribu¬ 
tion width, prestress and bulk and interface stiffnesses 
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FIG. 13: Color) Rupture front speed depends on local 
strength, (a) Rupture front speed vs position along the 
interface for a series of simulations in which the 
distribution of junction forces was varied systematically 
and two values of shear prestress were used. All of these 
events are of the fast-only type. The magenta (grey) 
and black lines have the shear stress values listed in the 
legend of (b), these also correspond to the magenta 
(grey) and black lines in Fig. 12c. The width a of the 
junction force distribution is increasing from the bottom 
to top line of each color; the values can be found in (b). 
For each value of junction force distribution width, the 
front propagation speed is higher when the prestress is 
higher, that is, each magenta (grey) line is higher than 
the corresponding black line, (b) Fast front speed data 
from (a) vs junction force distribution width. For this 
dataset we define, for each prestress series, the front 
speed at cr = 0 as the reference speed Uref- Then, for 
each simulation we calculate, for every position, 
{vc{x) — Vref{x))/vref{x). The dots are the mean of these 
values, the bars show the min/max values. The values 
on the cr/(/thres — /slip) axis range from 0 to 0.41 in 14 
equal steps of 0.029. (c) Slow front speed vs junction 
force distribution width. The region over which the 
slow front extends was determined from visual 
inspection of the underlying data, which is shown in 
Fig. 21. The prestress is tq = 0.3 in the triggering 
region and tq = 0.05 in the propagation region; the 
junction distribution widths can be found in (d). (d) 
The slow front speed data in (c) vs junction force 
distribution width. The same colors/grey levels are 
used in (c) and (d). Dots: mean values. Bars: min/max 
values. The data values on the cr/(/thres ~ /slip) axis 
range from 0.12 to 0.47 in 6 equal steps of 0.059. 


(Fig. 14c shows the same data on logarithmic axes). 
These blocks are all in the propagation region and away 
from the triggering zone and the leading edge. The front 
speed at each point along the interface is measured in 
the same way as in Fig. 3. The corresponding slip speed 
is the average block slip speed measured from the time 
the block, i, starts slipping to the time when the block 
5 mm away (block i -\- 2) starts slipping. We choose this 
time interval because it is during this period that the 
motion of block i most directly affects the rupture front 
propagation. 

Figure 14b shows the data in Fig. 14a with the fast 
slip speed in each simulation rescaled according to the 
parameters used (Fig. 14d shows the same data on loga¬ 
rithmic axes). We also include (in black) the slow front 
and slow slip data from [27, Fig. 3D]. Arguments for 
the scaling were presented in [27, SI Equations]; we re¬ 
peat only the conclusions here. The rescaled slip speed 

is Ugiip^ rescaled ~ '^slip rthres— tq ’ W^-CrC kij is the 

stiffness of the connection between a block and the in¬ 
terface, lo = 7 mm is the characteristic decay length of 
the shear stress field, which depends on the bulk to in¬ 
terfacial stiffness ratio k/ki^ Tthres = A^s/thres, and tq 
is the initial shear stress before the event. The normal 
force enters in the scaling indirectly because it modifies 
/thres- We expect deviations from the rescaling due to 
two sources that have not been included in the scaling 
equation. First, the effective force threshold on a 

block with a given junction force distribution is always 
less than Tthres? as discussed in Appendix B. The data 
where the junction force distribution was varied collapses 
better when replaces Tthres ia the scaling relation. 

However, because T^j^^g depends on the microscopic state 
of the interface, we choose to omit this correction in order 
to keep the factors in the scaling equation at the meso¬ 
scopic level. Second, the shape of the shear stress field 
also appears in the argument for the scaling, but not in 
our final scaling equation, which only includes the char¬ 
acteristic stress decay length Iq; we have kept its variation 
between simulations to a minimum by keeping the ratio 
k/ki constant. 

Figure 14b shows that the above rescaling allows all 
data from Fig. 14a to nicely collapse on a single straight 
line, going through the origin and having a slope very 
close to 1. This demonstrates that the front speed is 
directly proportional to the concurrent slip speed, with 
the proportionality coefficient being ——. Strikingly, 
the data for slow slip and slow fronts collapse on the very 
same line, as clearly seen in Fig. 14d. We return to this 
result in the discussion. Section VII. 


E. Front speed compares well to experiments 

So far in Section V we have sought to isolate the influ¬ 
ence on front speed of transients, prestress and interface 
strength from each other and from the event triggering. 
We now combine these results in a form suitable for com- 
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Slip speed, (m/ s) Rescaled slip speed, (m/ s) 

FIG. 14: Color) Front speed is proportional to slip 
speed, (a) Front speed vs slip speed from simulations 
with controlled initial states in the regime that 
produces fast fronts. Taking the simulation that plots 
at (0.23, 0.2) in Fig. 9 as an arbitrary reference, we 
have varied the junction force distribution 0 (/t) 
(magenta crosses), the prestress fo (green dots), the 
mass density p (blue squares), the bulk and interfacial 
stiffnesses (keeping k/ki constant) and p (cyan circles), 
and the normal force F/v and p (red diamonds), (b) 
The data in (a) with slip speed rescaled as 

'^siip, rescaled = '^siip rthre^-ro ’ markers near the 

origin show the data for slow slip and slow fronts from 
[27, Fig. 3D], (c) and (d) The data in (a) and (b), 
respectively, on logarithmic axes. 


parison with experiments. 

Figure 15 shows the rupture front speed vs prestress ra¬ 
tio r/p for several events where the shear prestress and 
the junction force distributions were varied systemati¬ 
cally. In order to both enable the tuning of prestress and 
force distributions and to mimic the experimental driv¬ 
ing conditions more closely, as in the full simulations, 
initial conditions were prepared in two steps. First, ho¬ 
mogeneous distributions and prestresses were assigned. 
Then, the trailing edge driving spring was allowed to 
move, eventually triggering the front. This raises the 
prestress most prominently near the trailing edge, but 
for the present parameters (which are found from com¬ 
parison with experiments) the ratio of shear stress profile 
decay length, /q, to system length is such that the shear 
stress is non-negligibly changed everywhere. For each 
event, rupture speed is measured at all points between 
X = 2.5 cm and x = 11 cm, i.e. avoiding the immediate 
vicinity of the system edges. 



Stress at front start, r/p 

FIG. 15: Golor) Rupture front speed vs prestress for 
several events where the prestress and the junction force 
distribution have been varied systematically. Different 
colors correspond to different widths of the junction 
force distribution: cr/(/thres ~ /slip) = 0.014 (red), 0.143 
(green), 0.286 (blue); the bell-shaped distributions 
defined in Fig. 20 were used. Different markers 
correspond to different initial levels of prestress (before 
additional loading was applied): r = —0.02 (t), 0.0 (•), 
0.06 (a), 0.13 (■), 0.2 (◄), 0.3 (★), 0.4 (►), 0.5 (♦). In 
terms of r/p = fsiip/fN + T(/thres//jv - /siip//jv), these 
prestresses are 0.165, 0.170, 0.184, 0.200, 0.216, 0.239, 
0.262 and 0.285, respectively. Edge effects, which can 
dominate visually, have been excluded by showing only 
the region between x = 2.5 cm and x = 11 cm. 
Prestress is measured for all blocks when the first block 
starts slipping (as defined from the number of pinned 
junctions). To get a sense for the range of values on the 
r/p axis, which differs from the ones in [16], consider 
that values much below /sUp/Zat = 0.17 lead to 
arresting fronts and that the maximum possible value of 
fibres/I n = 0.4 is reached only for narrow junction force 
distributions and near the loading point, so it is 
excluded by excluding the edges of the system. Recall 
from Fig. 9 that combinations of high prestress and 
high a are forbidden: this is the reason why there are 
fewer blue entries. For these simulations, a time step of 
At = 5 • 10“^ s was used to obtain smooth 
measurements at high front speeds. 


Figure 15 is analogous to Fig. 3 in Ben-David et al. [16], 
which shows the rupture front speed in selected points 
away from the sample edge vs prestress. Both figures 
show the same trend, namely that the front speed in¬ 
creases with increasing prestress, and that there exists 
a continuum of front speeds ranging from close to zero 
(slow) up to super-shear wave speeds. In addition, the 
numerical data indicates a systematic effect of varying 
the width of the junction force distribution, in a way 
which is completely consistent with the observations re¬ 
ported in section V G: larger widths yield lower local 
strength and thus promote faster fronts. 
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VI. HISTORY DEPENDENCE OF MESOSCOPIC 
STATIC FRICTION COEFFICIENT 

In the previous section we argued that the rupture 
front speed depends on the local strength, i.e. the 
effective static friction threshold on the block level. 
Here we investigate the extent to which variations in 
the mesoscale strength occur even when the individual 
micro-junctions attached to the block all have the same 
strength. This effect comes in addition to spatial hetero¬ 
geneities in stresses, due for instance to previous ruptures 
of the interface, which also cause spatial variations in the 
mesoscale strength. 

We showed in [47] (and recall in Appendix B) how 
the local strength depends on the distribution of forces 
among the micro junctions of a block. We also showed 
how, for simple velocity profiles, the distribution of forces 
at block arrest depends on the deceleration that led to 
this arrest. More fundamentally (in the model), the dis¬ 
tribution of forces at arrest depends on the slip profile 
x{t) of the block during the re-pinning of its junctions. 
As simplified velocity (and thus slip) profiles were studied 
in [47], we focus here on the more complicated sponta¬ 
neously occurring events in the full simulations. We are 
primarily interested in the effective static friction thresh¬ 
old = Tmax/p after an event, where Tmax is the largest 
tangential force that the block will bear before it starts 
to slip in the next event, can be calculated directly 
from the junction force distribution if we assume that the 
rupture is fast compared to the mean re-pinning time. As 
we will show, however, can also be estimated from 
the preceding slip dynamics. 

This section has four figures where the first three serve 
to introduce new concepts and definitions and the fourth 
is the main result of the section. Figure 16 defines the 
Gini coefficient, which we will use as a scalar, integrated 
quantifier of complicated block slip dynamics. Figure 17 
provides a link to [47, Section IV] in using the same con¬ 
stant deceleration slip profiles as were used there. These 
will also serve as the basis for our prediction of effective 
static friction as a function of Gini coefficient. Figure 18 
shows example slip data from full simulations and how 
we measure the Gini coefficient from these data. Finally, 
Fig. 19 shows that the Gini coefficient is a good predictor 
of effective static friction. 


A. Gini coefficient: definition and measurement 

A feature of the microscopic junction law that we 
also studied in [47] is how the junction force distribu¬ 
tion 0 (/t) of a block evolves with the slip history of the 
block. As the block starts moving, slips and comes to 
rest, junctions are continually breaking and reforming, 
and the complete distribution 0 (/t) depends on the full 
slip velocity history v{t). In practice, however, it is the 
very end of the preceding event that has the largest in¬ 
fluence on 0, and this allows us to predict the strength 


100 % 



FIG. 16: The areas that define the Gini coefficient, a 
measure of inequality. 


of (j) after an event (which is the effective static friction 
threshold for the next event) from a characterization of 
the block slip motion. 

First, we pick a time t' between two full sliding events 
at which we wish to predict the strength of (j). The blocks 
should be at rest at this point in time. Second, we iden¬ 
tify the segment [to,ti] of the slip history that needs to 
be considered to predict the strength of (j) at the chosen 
time. Following each block’s motion x{t) backwards from 
t', we define to as the time when the block was one junc¬ 
tion breaking length 5m = fthres/hj away from where it 
is at t'; x(t') — x(to) = 5m- The slip history before to can 
be neglected because all the junction attachment points 
that existed before to have been broken and renewed dur¬ 
ing t G [to,t']. Further, we define ti = to + ti? + 6tR, so 
that [to,ti] includes most of the junction renewal occur¬ 
ring after to- We require ti < t', but as long as the blocks 
remain at rest, any ti that satisfies this constraint can be 
chosen. It is possible to create pathological cases where 
[to,ti] fails to include the bulk of re-pinning events that 
set up 0 (/t) at t', but in practice, in all our simulations 
this is the most important segment of the slip history. 
Third, we characterize v{t) for t G [to,ti] and compare 
it to velocity profiles that can be treated analytically or 
with simple numerical integration. 

The feature of v{t) that has the most prominent im¬ 
pact on (j) is how evenly the slip motion x{t) is distributed 
in time. If the block comes fully or nearly fully to rest 
during a time interval much shorter than tR, most of 
the junctions will be renewed while the block is at rest; 
they will therefore have approximately the same stretch¬ 
ing and 0 (/t) will be narrow. If the block moves more 
uniformly as it comes to rest, the stretching of the junc¬ 
tions will distribute more uniformly. This simple argu¬ 
ment assumes that the rate of junction renewal is con¬ 
stant, a simplification that is a reasonable approximation 
when the block comes to rest after significant slip motion 

[47]. 

A robust measure of inequality is the Gini coefficient, 
which we introduce in Fig. 16. It is commonly used to 
characterize the inequality of the income or wealth dis¬ 
tribution in a population [60, 61], but can be applied to 
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our case without modification. The Gini coefficient of a 
set is defined from the areas in Fig. 16 as ^ = Aj 
(which equals 2A, since A B = 1/2). For any set of 
elements, each having a share of a some quantity (e.g. 
wealth), the line separating A and B is constructed as 
follows. Elements are first ordered by the size of their 
share and placed along the abscissa. The ordinate for 
a given share then corresponds to the sum of all shares 
that are smaller than the chosen share. The line thus 
defines the cumulative distribution of shares. B is the 
area under the cumulative and A is the area between the 
cumulative and the line of equality. The line of equality 
is the cumulative for a set where each element has the 
same share. Thus, if each element has the same share, 
0 = 0, while in the other extreme where one element has 
everything and the other elements have zero share, 0 = 1. 
To measure the Gini coefficient of the block slip during 
t G [to, ti] we create a set consisting of Atg = 10 elements. 
We define Ateiem = (ti — to)/Ne and assign Xghare to 
each element so that the elements tile the time interval. 
That is, the share of the first element is the block slip 
that occurred between to and to + Ateiem, the share of 
the second element is the block slip that occurred during 
[to + Ateiem, to + 2Ateiem], and SO on. Equivaleutly, each 
element’s share is the average slip velocity during the ele¬ 
ment’s time interval (net slip equals average slip velocity 
times the length of the time interval, but the common 
normalization factor Ateiem can be ignored). 

Eigure 17a shows block slip profiles x{t) for constant 
deceleration slip. These slip profiles have the advantage 
that they are completely characterized by the decelera¬ 
tion amplitude; they help build intuition before we con¬ 
sider the full simulation data, where none of the basic 
kinematic quantities (slip, velocity, acceleration) are even 
approximately constant. In Eig. 17a the time interval has 
been shifted so that to = 0. With the parameters in Ta¬ 
ble I, (ti — to)/tR = 1.3. Eigure 17b shows the effective 
static friction, i.e. the strength of (l){fT), resulting from 
the slip profiles in Eig. 17a. We showed in [47, Eig. 9] 
that for (i) initial velocities large enough that the blocks 
move longer than the junction breaking length be¬ 
fore coming to rest, and (ii) reasonable initial junction 
distributions, the effective static friction for constant de¬ 
celeration slip depends only on the deceleration value and 
the junction law parameters; here we use the same initial 
conditions as in [47, Eig. 9] and the junction law param¬ 
eters in Table I. Eigure 17c shows the Gini coefficient of 
the slip profiles in Eig. 17a as a function of their constant 
deceleration parameter (Gini coefficients were calculated 
using [72]). The shape of the curve can be intuitively 
understood. Eor deceleration close to zero, the slip mo¬ 
tion from to to ti occurs with low and nearly constant 
speed, so the cumulative is close to the line of equality, 
and 0^0. Eor large deceleration values, the block comes 
to rest sharply, [to,ti] includes a long period where the 
block is at rest, and the cumulative is more peaked. In 
the limit of very large deceleration values, 0^1. We 
can combine Eig. 17b and Eig. 17c to obtain a prediction 




EIG. 17: Golor) (a) Slip profiles x{t) for constant 
deceleration slip. The markers extend over the interval 
[to,ti]. (b) Rescaled effective static friction for the 
slip profiles in (a) and the junction parameters in 
Table I. (c) Gini coefficients of the slip profiles in (a). 
The line extends to the origin. The markers in (b) and 
(c) show the data from the slip profiles in (a) of the 
corresponding color/grayscale; the drawn lines are 
based on a denser set of slip profiles that were omitted 
from panel (a) for clarity. Combining (b) and (c) gives 
the drawn blue line with circular markers in Eig. 19. 


for vs 0. This gives the full drawn lines in Eig. 19. 

Eigure 18 shows example slip profiles from full sliding 
events in full simulations. Two events are shown, with 
the data in c being a detailed view of the data in a and 
the data in d being a detailed view of the data in b. 


B. Gini coefficient predicts effective static friction 

We are now ready to test how well the Gini coefficient, 
an integrated property of the block slip history, predicts 
the effective static friction threshold. The link between 
the two is the distribution of junction forces that exist 
after a sliding event: the Gini coefficient is an averaged 
measure of how uniform or non-uniform the distribution 
becomes; the effective static friction threshold is essen¬ 
tially a property of the distribution (it also depends on 
the junction evolution law, but this remains fixed within 
a simulation). 

Eor the limiting cases of a block coming to rest very 
abruptly or very gradually, the Gini coefficient and the 
effective static friction threshold can both be predicted 
exactly from simple arguments. The abrupt stopping 
case is the limit where the block comes to rest from a 









19 


a 





Time (ms) 


Time (ms) 


FIG. 18: Color) Example slip profiles from full 
simulations, (a) and (b) Slip profiles for full sliding 
events. Only the slip occurring between to and ti 
(defined in the text) are used to calculate the Gini 
coefficient. Effective static friction is measured from the 
junction force distribution (/)(/t) at t'. The data in (a) 
are for the block at x = 4 cm in Event I, while the data 
in (b) are from a simulation with /siipZ/w = 0.05. (c) 
and (d) Detailed views of the slip profiles in (a) and (b), 
respectively, for the time interval between to and ti. 
Black line without markers: x{t). Colored line with 
markers: a moving average filter of width 1 ms was 
applied to x{t) and then Ne = 10 equally sized intervals 
were defined (markers indicate intervals’ boundaries). 
The smoothing helps avoid negative increments that 
while possible to include are not part of the basic 
formulation of the Gini coefficient. 


is large compared to the junction breaking length 
the motion is always far from the abrupt stopping limit. 
Consequently, to span out the range of Gini coefficients 
from 0 to 1 we have performed simulations where we in¬ 
directly varied the relative amplitude of slow slip and 
by varying /thres and /sUp. Eigure 19 shows results from 
this series of full simulations. As expected, the variation 
between events within one simulation is small compared 
to the variation between simulations. In all cases, the 
Gini coefficient is a good predictor of effective static fric¬ 
tion. In particular, we find near equality between Q and 
for Q > 0.5. 

The agreement between the line corresponding to con¬ 
stant deceleration slip and the simulation data is rather 
good, but far from perfect, as expected given the approxi¬ 
mations inherent in describing the complicated block slip 
motion arising in the full simulations with a single num¬ 
ber. Important differences in the underlying motion are: 
(i) In the full simulations the duration of sliding events 
is on the order of tR, see Eig. I8a and b. To reach steady 
state junction distributions, as was assumed in [47], a 
longer sliding period is necessary, (ii) In the full simula¬ 
tions, the oscillations in the block velocity set up junction 
force distributions that are non-smooth. This is why the 
simulation data is not bounded by the same lower limit 
as the constant deceleration data. 


VII. DISCUSSION 

We begin this section with an overview of the robust¬ 
ness of the observed front dynamics to variations in model 
parameters. We then revisit the results and provide in¬ 
terpretations and comparisons to earlier work in the lit¬ 
erature. 


A. Robustness 


speed high compared to s^/tR within a distance short 
compared to Sm and a time short compared to Ir. In 
this case the entire motion between to and ti will be 
assigned to a single element and the Gini coefficient 
will he Q = I. The distribution of junction forces will 
be a (^-function as all junctions reform after the block 
has stopped moving, and so Tmax/p = /thres//iv and 
the rescaled effective static friction threshold is /i®® = 

(Tmax/p-/siip//iv)/(/thres//jv-/siip//jv) = 1- The grad- 
ual stopping case is the limit where the block spends 
a long time at velocities small compared to Sm/ti?. In 
this case the velocity is nearly constant from to to ti 
and the Gini coefficient will be ^ = 0. A uniform junc¬ 
tion force distribution is set up (see [47, Section III A] 
for a detailed account). Erom equation (B2) we find 

Pf = (/thres - /slip)/(2/thres)- 

As a slipping block slows down the slow slip mech¬ 
anism becomes important. If the amount of slow slip 


To elucidate the robustness of the system dynamics to 
variations in the model parameters we have performed 
simulations where one or two parameters at a time were 
varied from their values in Table I. We find that most of 
the parameters are in a range where the results remain 
qualitatively the same under halving and doubling of the 
parameter value. Others, notably the bulk and interface 
compliances, are in a range where these changes lead to 
significant changes in the system dynamics. 

As mentioned in the model description in Section II, 
T(tR) is a simplified way of modeling the distribution 
of times after which slipping micro-junctions relax. It is 
therefore reassuring that the front dynamics are not sen¬ 
sitive to the details of T(tR)\ (i) We observe fast-slow- 
fast fronts in a simulation with an exponential T(tR). (ii) 
In [27] we varied the average, {Ir)^ systematically and 
observed that this changes the speed, but not the na¬ 
ture of slow slip and slow fronts, (hi) We have performed 
full simulations where we halved and doubled the width. 
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FIG. 19: Color) Rescaled effective static friction after 
full sliding events vs Gini coefficient, an integrated 
quantifier of block slip motion. All data from full 
simulations, with variation in junction law parameters 
between simulations as shown in the legend. See Fig. 17 
for an explanation of the Gini coefficient and the 
definition of the theoretical lines, and Fig. 18 for sample 
block slip histories. For each simulation we show data 
from two arbitrary full sliding events (solid and open 
symbols). 


StR^ from its value in Table L Although the details of 
the fronts change, we still observe precursors, fast-only 
events, fast-slow-fast events and intermediate partial slip 
events. 

The number of junctions per block, Ng, could be ad¬ 
justed to match the area density of junctions in a given 
experiment. However, for ^ 1 we expect to recover 
the results from [47], where instead of tracking individ¬ 
ual junctions we described the state evolution in terms 
of continuum distributions of pinned and slipping junc¬ 
tions. With Ns = 100, the present simulations and those 
in [27] appear to be within this continuum limit: when 
we compare to full simulations with Ng = 50, 200 and 
500 we again observe precursors, fast-only events, fast- 
slow-fast events and intermediate partial slip events, and 
although due to the randomness in the junction dynam¬ 
ics the events are not exactly the same between these 
simulations, they are very similar. 

There are four stiffness parameters in the model: the 
driving spring modulus K] the elastic foundation modu¬ 
lus kf] Young’s modulus from which the bulk spring 
stiffness k follows; and the interface spring stiffness kij. 
We have performed simulations where the driving spring 
modulus K was increased twofold, tenfold and hundred¬ 


fold from its value in Table L This changes the relatively 
regular stick-slip seen in the loading curve (Fig. Id) to 
a more chaotic stick-slip pattern, and as predicted in 
Section IV B it causes the amplitude of the drops in Ft 
during partial slip events to increase. We still observe 
both fast-only and fast-slow-fast fronts. For the elastic 
foundation modulus, kf^ doubling and halving its value 
in full simulations still allows the full range of system 
dynamics. For kf k, the normal stress profile p{x) is 
essentially fiat, as significant variations in p would then 
require large vertical deformations of the slider, and the 
experimental stress profiles of [6, 16] will no longer be 
reproduced. 

For the bulk and interface stiffnesses k and kij we have 
shown in Fig. 14 and in [27] that when their ratio remains 
constant, changes to the value of these parameters affect 
the slip and front speed, but dynamics remain qualita¬ 
tively the same. When their ratio changes, however, qual¬ 
itative changes to the system dynamics occur. Doubling 
k or halving kij increases the loading region where the 
interface stress is significantly influenced by the increas¬ 
ing force in the driving spring between events. This leads 
to multiple front tips at different locations propagating 
simultaneously instead of one front propagating from the 
trailing to the leading edge, so that a single front speed 
cannot readily be defined. Halving k or doubling kij^ on 
the other hand, reduces the extent of the loading region, 
so fronts can still be defined as before, but these settings 
favor fast front propagation, so very few fast-slow-fast 
fronts were observed. Of course, the slow slip mechanism 
is still active, but as discussed before it is masked by the 
fast slip and fast front propagation. 

Variations in the force levels in the junction evolution 
law, /thres and /slip, are included in the scaling results 
that are explained in Appendix B. We also fully expect 
that parameter values can be found where full simula¬ 
tions produce precursors, fast-only events, fast-slow-fast 
events and intermediate partial slip events with larger 
values of /thres and /sUp than we use. Because this change 
modifies both the triggering and propagation of events, 
other parameters would have to be changed simultane¬ 
ously. 


B. Interpretations and comparisons 

In Section IV B we identified a signature of slow front 
propagation in the macroscopic loading curve. If this 
signature turns out to also be present in the experiments 
where slow fronts are known to occur, it could be a use¬ 
ful first indicator of slow fronts in experiments where the 
nature of rupture front propagation is unknown. The 
underlying insight is that the reduction in Ft during 
an event has two distinct parts, both of which are due 
to the motion of the trailing edge of the slider: first. 
Ft drops during front propagation, as the slider deforms 
when some parts of the interface are slipping and others 
are stuck; then. Ft continues to drop during full sliding. 
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If the reduction in Ft during front propagation is suffi¬ 
ciently large to be reliably measured, a fast-slow transi¬ 
tion will show up in as a reduction in slope. If the 

front transitions back to fast propagation or reaches the 
leading edge, a change back to large negative slope will 
occur in Frit). 

In Section IV C we demonstrated by changing the slider 
mass density p and hence the bulk wave speeds that the 
fast slip and fast fronts in the model are of inertial origin. 
This result, which was to be expected, further supports 
our claim from [27] that fast slip and slow slip are due 
to distinct mechanisms. It follows that the ratio of their 
speeds can vary if system parameters are changed. For 
instance, in the model: fast slip speed is independent of 
and slow slip speed is inversely proportional to the time 
scale Ir which determines interface healing after rupture 
[27, equation SI]; fast slip speed is inversely proportional 
the square root of slider mass density and slow slip speed 
is independent of slider mass density. In contrast, reduc¬ 
ing (increasing) the interface stiffness increases (reduces) 
both slow and fast slip speeds and thus does not affect 
their ratio. Because front speed is proportional to slip 
speed across the velocity scales that we have observed, 
see Fig. 14, it follows that large variations in the ratio of 
slow to fast front speeds is also possible between systems. 
This paradoxically suggests that the usual definition of 
slow fronts in terms of their velocity (one to several or¬ 
ders of magnitude slower than the Rayleigh wave speed), 
although natural in the system in which they were dis¬ 
covered [3], is not satisfying. This definition does not 
capture the physical characteristics of slow fronts - dy¬ 
namic fronts propagating due to the intrinsic relaxation 
dynamics of the frictional interface after arrest - or dif¬ 
ferentiate them from other front types, and may lead to 
misinterpretations in systems where the velocity range of 
slow fronts significantly overlap the velocity range of fast 
fronts (of the order of sound speed) or quasi-static fronts 
(speed proportional to the rate of external loading). 

In Section IV F we first repeated our finding from [27] 
that the type of front observed (arresting, fast-slow-fast, 
or fast-only) depends systematically on both the strength 
and the prestress of the interface, with weaker interfaces 
and higher prestresses favoring faster front propagation. 
We then showed that although the front type phase di¬ 
agram in Fig. 9 is not a local measure, because front 
speed is transient and also depends on the interface state 
behind the front tip, the diagram can still give useful 
predictions of the way front type changes with sudden 
changes in interface conditions. This could help in the 
interpretation of experiments on interfaces with hetero¬ 
geneous friction properties like [73], where material het¬ 
erogeneities in the form of rock pebbles embedded at a 
gel-glass interface were observed to cause the rupture 
front speed to sometimes increase, sometimes decrease. 

In Section V B we showed that, as expected from exper¬ 
iments and other numerical work, front speed increases 
with increasing shear to normal prestress ratio r/p. How¬ 
ever, knowing the local r/p alone is insufficient informa¬ 


tion to predict the front speed. Indeed, in Section V A we 
showed that in the model and for experimentally realistic 
system sizes, front speed is transient at all points and as 
such the selection of its value is both history-dependent 
and non-local. In Section V C we showed that both fast 
and slow front speeds also depend on the local strength 
of the interface, which is generally out of experimental 
reach, especially when the interfacial junctions bear some 
randomness, e.g. in their size or stiffness. 

In [27] we demonstrated that the slow front speed in 
the model is proportional to the slow slip speed and 
worked out the constant of proportionality. Bar Sinai 
et al. [29, equation (9)] (see also [28, equation (9)]) found 
a similar result in an analytical model with a rate-and- 
state-type friction law. In Section VD we found that 
the same constant of proportionality relates fast slip to 
fast front speed. The proportionality follows from relat¬ 
ing how far the blocks behind the front tip have to move 
to advance the rupture, to how quickly they move (the 
slip speed). For instance, less motion is required, lead¬ 
ing to faster fronts, when the interface is stiffer or more 
highly prestressed, or when the stress field in the bulk 
has a longer decay length. Taken together, these results 
demonstrate that the origin of the front speed, be it slow 
or fast, lies in the relevant slip speed during front prop¬ 
agation. The selection of a fast or a slow front regime is 
thus the direct consequence of the selection of a fast (in¬ 
ertial) or a slow (driven by the internal relaxation of an 
arrested interface) slip regime. We therefore claim that 
understanding front type selection and the transition be¬ 
tween front types essentially reduces to understanding 
the slip motion and its time evolution during an event. 
Note, however, that even when a mechanism for slow slip 
is present in a system, slow fronts will not appear if the 
fast fronts always traverse the entire interface: being fast, 
they will mask any slow propagation. 

In Bar Sinai et al. [29], the slow slip speed was se¬ 
lected as the speed for which the steady-state friction 
law had a minimum. However, the existence of such a 
minimum is not a necessary condition for slow fronts to 
appear in a system. As a matter of fact, in our model, 
and for the parameters used, the steady state friction law 
is purely decreasing. The apparent contradiction disap¬ 
pears as soon as one remembers that in our model, the 
slow slip speed arises from a completely different phys¬ 
ical mechanism, namely the time distributed relaxation 
of the interface after its rupture and arrest. 

For the scaling of the shear to normal stress ra¬ 
tio we have used tq = (tq/p - fs\ip/I n)/{ fthres/I n - 
/sUp/Zat) throughout. It would have been consistent with 
the observed dependence of front speed on both local 
stress state and local strength to replace fthres/fN with 

Tmax/p = to obtain = (ro/p-Zslip/Zw)/(Tmax/p- 

/sUp/Zat), and use this in Fig. 12b. Nevertheless, we have 
chosen to stick with a scaling relation that depends only 
on input model parameters and not on the emerging ef¬ 
fective strength of the interface, because quantitative val¬ 
idation of the more advanced scaling is not directly pos- 
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sible due to the transient nature of the fronts. In an ap¬ 
plication where fronts were known to have converged, it 
would be interesting to try the scaling when plotting 
the fronts’ stationary speed against prestress. Similarly, 
replacing Tthres with Tmax in the scaling relation between 
front speed and slip speed in Fig. 14 would improve the 
scaling of the data where 0 (/t) was varied, at the cost 
of introducing emerging properties into the scaling. 

Ben-David et al. [16, Fig. 3] reported a relationship 
between front speed and shear to normal prestress ra¬ 
tio. In Section V E we showed that although our results 
are quantitatively different (in particular, the range of 
tq/p values that are stable in the model is about half 
of that in the experiments), qualitatively we capture the 
experimental observation well. Further, our results for 
front speed transients (Section V A) and dependence on 
local strength (Section V C) are possible explanations 
for the relatively large range of front speed values that 
were observed for each value of r/p in [16]. Namely, 
we expect variation in local strength to be present in 
the experiments, probably to a larger degree than in the 
simulations, as we have explicitly avoided spatial hetero¬ 
geneities in the friction parameters. It is also possible 
that the experimental fronts share with our simulations 
the property that over the sample length studied, the 
front speed is in the transient regime: the fronts in e.g. 
Fig. 2 of [16] actually change their speed throughout the 
interface, but whether this is due purely to changing in¬ 
terface conditions or also to a “true” transient behavior 
is difficult to assess. 

We argued in the model description (Section II) that 
a time-dependent rule for the reformation of junctions 
has experimental justification. A consequence of this mi¬ 
croscopic friction law is that the local microscopic state, 
specifically the distribution of forces in the junctions, de¬ 
pends on the block slip dynamics of the preceding event 
[47]. The microscopic junction state in turn determines 
the effective static friction, which affects the front propa¬ 
gation. In Section VI we showed that the Gini coefficient, 
an integrated estimator of the non-uniformity of the block 
slip history, predicts the effective local static friction even 
in full simulations were the block slip dynamics is highly 
complicated. We were thus able to link two very differ¬ 
ent aspects of the system behavior. Further, while the 
distribution of junction forces and the resulting effective 
friction are difficult to measure experimentally, the local 
slip history can more readily be measured, at least on the 
side of the slider [see e.g. 22, 59] or in setups which track 
the motion of patterns or markers at the interface [see 
e.g. 4, 5]. We recognize that there are challenges that 
need to be overcome to adjust our protocol for measur¬ 
ing the Gini coefficient and apply it to experiments or 
real systems like seismic faults. However, if the length 
scale for junction breaking and the time scale for junc¬ 
tion reformation can be determined, we believe that the 
underlying idea, namely to characterize the crucial part 
of the slip motion with a robust measure of inequality 
and to use it to predict the (still unknown) strength of 


the interface, is applicable. 

We have studied the onset of sliding friction in a model 
that couples a microscopic friction law involving the state 
of a large number of individual junctions to a 2D elastic 
solver. In the range of possible combinations of friction 
and bulk models outlined in the introduction this is, we 
believe, an example of bringing together a moderately 
complicated friction law with a moderately complicated 
bulk law. For example, the friction law excludes both ag¬ 
ing and a distribution of junction strengths, and the bulk 
law, while 2-dimensional, excludes among other things 
plasticity and melting, and is solving for the deformation 
of the slider only. Nevertheless, the simulations repro¬ 
duce many of the features of the spatio-temporal dynam¬ 
ics in experiments to which the model was attuned, and 
we were able to use the more complete information of the 
system and the fine control of initial conditions available 
in simulations to better understand how some of the cen¬ 
tral material properties and system state characteristics 
modify the dynamics. We believe important avenues for 
future research to be the inclusion of the bulk deforma¬ 
tions of the substrate, which brings with it the compli¬ 
cation of defining frictional properties at a non-planar 
interface, and the parametrization from more fundamen¬ 
tal models of the individual junction evolution laws. 
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Appendix A: Stencil for determining rupture front 
speeds 

In this appendix we give the details of how we measure 
the front propagation speed Uc- We define it as 


where Xtip(t) is the location of the front tip. To use this 
expression in practice, the position of the front tip in time 
must be found and a suitable stencil used to take the 
derivative numerically. We use the rupture criterion on 
the block level, that less than 30% of a block’s junctions 
remain in the pinned state, to find rupture times trup,i 
for every block that partakes in an event. We also define 
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the front tip at any time to be located at the last block 
to reach this sliding criterion. The simplest formula for 
the local rupture speed is then to take the average speed 
of the front on its way between two blocks, that is 


Xi^l - Xi 


t 


rup,i+l 


- tr 


(A2) 


However, this is very sensitive to the discreteness and 
non-smoothness of the front and gives very large fluctua¬ 
tions along the interface. Smoother results are obtained 
by averaging over a larger area. We have used 


^i+vwidth vwidth / a 

= - - — -, (A3) 

^rupji+vwidth ^rup,i—vwidth 

with vwidth a number of blocks. We have found that 
vwidth = 2, which gives a 5-point wide stencil, strikes a 
good balance between smoothness and spatial resolution 
for the parameters we employ. Of course, for higher 
Nz, vwidth can be increased proportionally. 

In principle, another option for finding a larger kernel 
is to use standard stencils for numerical derivatives of 
first order. Since the raw data are the progressive trup,i 
over an approximately constantly spaced spatial grid, the 
inverse rupture speed can be found from the standard 5- 
point stencil as 


^rup,i-|-2 T 8^rup,i-|-l —1 T ^rup,i —2 

12Ax ' 


(A4) 


Although this is a kernel with five nodes, the weights 
are such that the nodes near the centre dominate, and 
its fluctuations are comparable to those for vwidth = 1 
in Equation (A3) and larger than for vwidth = 2. The 
standard 7- and 9-point stencils are also dominated by 
the centre nodes, giving poor results. 

When a front nucleates, either as a first front nucleat¬ 
ing in a fully pinned interface, or as a secondary front nu¬ 
cleating ahead of the main front due to inhomogeneities 
in the stress state, it travels both left and right. This 
leads to points on both sides of the nucleation site hav¬ 
ing nearly the same rupture time trup, and to the spuri¬ 
ous conclusion by a naive measurement scheme that the 
front speed was suddenly very high. To simplify auto¬ 
matic front speed measurement we have chosen to ignore 


any points for which the stencil in use extends across a 
region where the front is left-travelling. This is the rea¬ 
son for the front speed being undefined for some positions 
in Fig. 3c. 


Appendix B: The dependence of effective local static 
friction on the width of ^(/t), and the analytical 
result that underlies the scaling of and a axes 


In this appendix we show that the scaling that we use 
here and in [27] for the effective static friction thresh- 

old, jlf = {Tma^/p - /slip//jv)/(/thres//iV “ /slip//iv), 
and for the width of the junction force distribution, 
a = cr/(/thres ~ /slip), cau be derived analytically for 
uniform junction force distributions <p and holds also for 
other shapes of 0. The arguments were omitted from [27] 
for brevity. We include a figure that demonstrates the 
scaling and repeats the relationship between the shape 
and width of <p and the effective static friction threshold. 

The main result in this appendix is equation (B2). We 
will derive it using results from Thpgersen et al. [47], and 
some simple bookkeeping is required to set up the corre¬ 
spondence between our choice of parameters in that work 
and here. Namely, in the present paper and in [27] we 
express the state of the pinned junctions by the distribu¬ 
tion (j) of forces /t in the junctions. In [47] we expressed 
this state in terms of a length scale, the distribution S 
of junction stretching lengths 5. In the linearly elastic 
regime that we consider in all three works, these are re¬ 
lated by 0 = /cjun'S', where = kij is the single junction 
stiffness. 

In Appendix F of [47] we calculated the effective static 
friction for S uniform and with support of extent Sgupp- 
The assumptions were (i) linear springs with force 
at stretching 5, and (ii) breaking of the complete set of 
junctions attached to the block in a time interval small 
compared to the average time spent by a junction in the 
slipping state, so that the first contacts that break stay 
broken until is reached. Defining notation consistent 
with the present paper, we write equations (FI) and (F2) 
in [47] as 


Taax 


^2 (25thres '^supp) 

-Sslip -Sthres“l“'Ssupp ^ i 

g Jslip I 


2ss 


(s?hres 


'^slip 


, '^thres '^supp ^ '^slip 
, '^thres '^supp ^ '^slip, 


(Bl) 


where Sthres = = /thres/^jun Is the strctchiug at the junction breaking threshold and SgUp = /siip/^jun is the 

stretching corresponding to the sliding state friction. 

To rewrite equation (Bl) in terms of forces we define /supp = ^jun^supp- Subtracting /gUp and dividing by (/thres ~ 
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/slip) on both sides, and recalling that Ns/n = we arrive at 


-eff _ '^mayi/P /sUp/Zat 

/thres /In - /slip /In 


1 _ fs 

2 /thres 

1 res /slip 

2 /supp 


— /slip ’ /thres /supp ^ /slip 
5 /thres /supp ^ /slip- 


r 


(B2) 


On this form it is apparent that scaling by using the 
expression on the left hand side and scaling the width 
by using ^-f y uniform distribution with support 

of extent /supp has standard deviation a = /supp/(2\/3)) 
provides a data collapse under changes to /cjun, /slip and 

/thres- 

Figure 20a shows unsealed data for uniform (j). To 
show that the scaling in equation (B2) is not limited to 
this particular shape of 0, we include data for bell-shaped 
(j) in Fig. 20b. When we apply the scaling, in Fig. 20c, 
each distribution shape gets a data collapse. The master 
curves for the two distribution shapes are slightly differ¬ 
ent. 

A point which we have discussed in both of [27, 47] 
and which we use in Sections IV E and V C is that wider 
(/) correspond to lower effective local static friction. This 
can be seen directly in Fig. 20. Finally, we note that the 
scaling of is the same as the scaling of prestress f, a 
scaling that was also used in [9, 10, 74, 75]. 


Appendix C: Initialisation, boundary and driving 
conditions 


This appendix supplements the model description 
found in the main text with detailed information on 
how we initialize the system and apply the boundary 
conditions. The simulations we have performed can be 
grouped in three categories, (i) Full simulations. In these 
the sample is initially unstressed; then the normal force 
is applied; then the shear force is applied through the 
driving spring, and the spring drives the slider through 
tens of events: precursors, full sliding events, and pos¬ 
sibly partial slip events between the full sliding events, 
(ii) Simulations that restart at a point within a full sim¬ 
ulation. In these we make well-defined changes to iden¬ 
tify or highlight a particular mechanism, for example the 
slow slip mechanism, (iii) Smaller simulations where the 
shear and normal prestresses are controlled and a single 
event is studied. In these we have done triggering both 
through (iii)a breaking the junctions in a predefined re¬ 
gion, and (iii)b motion of the driving spring as in the full 
simulations. In this appendix we explain how each type 
of simulation was set up, and give parameter values in 
Table I. 
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FIG. 20: Color) Scaling of effective static friction vs 
junction force distribution width, (a) Unsealed 
analytical results for uniform junction distributions. 
Changes made to /cjun, /slip and /thres independently. 
Markers indicate /supp = /thres - /slip- (b) Unsealed 
numerical results for junction distributions that are 
bell-shaped polynomials with roots at ±a and the 
functional form 

HO = 5/(4a)(l +3|C/a|)(l - e [-a,a]. 

Changes made to /sUp and /thres independently, (c) 
Data collapse using exact scaling. 


1. Full simulations 


In full simulations, the slider is initialized with full 
normal load and no tangential load Ft by gradually 
applying Fn without allowing springs to break, a techni¬ 
cality required because the normal forces on the springs, 
fNiji start at zero and therefore springs, if allowed to, 
would break under any stretching. We distribute the 
load Fn uniformly on the top blocks; apart from this we 
use the same non-frictional boundary conditions as in [9]. 
The bottom blocks interact with a vertically elastic foun¬ 
dation that produces a normal force —kfZi when Zi < 0, 
zero when Zi > 0; Zi is the vertical displacement of in¬ 
terface block i; the value of kf is given in Table I. Zero 
force boundary conditions are applied to the right and 
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left edges (the top and bottom forces still apply to the 
corner blocks). The unique equilibrium is found through 
damped relaxation of typical duration 10 ms. After re¬ 
laxation, we check that no spring is stretched beyond its 
strength and introduce the driving spring starting from 
zero applied driving force Ft- Then Ft, which acts on 
the block on the left side of the slider situated at height h 
above the interface, through the driving spring, increases 
as the driving point moves to the right with speed V. 
Full simulations were used for Fig. 1, 2, 3, 4, 5, 18 and 
19. 


2. Restarted simulations 

It can be difficult to analyze events from the full simu¬ 
lations, because the initial conditions for each event arise 
dynamically from the preceding events. To make direct 
comparisons between events possible, restarted simula¬ 
tions begin from an interesting state arising in a full sim¬ 
ulation. The state consists of the instantaneous positions 
and velocities of all the blocks and the driving stage, and 
the pinned/slipping and lifetime information of all the 
junctions. We modify one or a few model parameters 
or state properties, rerun the simulation, and compare 
the output of the restarted simulation with the original, 
or with another restarted simulation. Restarted simula¬ 
tions were used for Fig. 6 (the mass density was varied), 7 
(junction evolution law was varied) and 8 (junction force 
distribution, a part of the state, was varied). 

3. Single event simulations 

Single event simulations, like restarted simulations, 
were performed to elucidate the qualitative and quan¬ 
titative importance of individual parameters or interface 
states by varying them individually, keeping the rest of 
the simulation setup unchanged. This allows us to go be¬ 
yond the insights from the full simulations. To simplify 
analysis, these systematic studies were done with differ¬ 
ent normal forces and different initialization from the full 
and restarted simulations. The normal force boundary 
conditions on the top and bottom were exchanged: this 
simplifies the analysis by setting a constant normal force 
Pi = Fn/Nx on all blocks i at the interface. To main¬ 
tain stability against global rotation, the top blocks in¬ 
teracted with an elastic ceiling with the same properties 
as the elastic foundation used in the full simulations. 

To obtain an initial state with a prescribed interfacial 
shear stress profile we turned the interface springs off 
during the initialisation. In their place we added to each 
bottom block the force corresponding to the shear stress 
to be prescribed. We also introduced the driving spring, 
but let R = 0. During relaxation, the sample moved 
along the x-axis until the force in the driving spring bal¬ 
anced the net force from the interfacial shear stress. To 
get rid of oscillations more efficiently we added damping 


forces —a{vi) on the blocks’ motion. After relaxation, 
the extra forces and the extra damping were turned off 
and the interfacial springs were introduced, with their 
attachment points Xij chosen such that the net force on 
each block was unchanged and the desired distribution 
of spring forces, 0(/t), appeared. We then waited a few 
timesteps to ensure that the transition from pre- to post¬ 
relaxation involved no force discontinuities. 

For type (iii)a single event simulations, instead of driv¬ 
ing the system with R / 0 until rupture was triggered, 
we started fronts by simultaneously de-pinning all junc¬ 
tions for all blocks to the left of ^trigger- The shear stress 
in the triggering region has a strong influence on the 
rupture fronts, and in order to simplify the comparison 
of results between simulations we used a constant value 
^trigger = 0.3. The triggering and propagation regions 
and the initial stress configuration for this type of sim¬ 
ulations are shown in Fig. 13c. Type (iii)a simulations 
were used for Fig. 9, 10, 11, 12, 13, 14 and 21. 

For type (iii)b single event simulations we matched 
the experimental event triggering conditions more closely. 
Again, homogeneous stress and junction force distribu¬ 
tion states were prepared as described above. Then, an 
event was started by moving the driving spring as in the 
full simulations. This modifies the shear stress state be¬ 
fore the event starts, particularly near the loading point. 
As in our other simulations, the width of the junction 
force distribution <p affects the effective static friction 
threshold and by that the energy released as the event 
starts. To isolate the effect of varying parameters in the 
propagation region, we kept the width of the junction 
force distribution constant at aioad/ (/thres~/siip) = 0.143 
to the left of xioad = 22.5 mm. Type (iii)b simulations 
were used for Fig. 11 and Fig. 15. 

For Fig. 11, which includes both type (iii)a and (iii)b 
simulations, we used the procedures described above, ex¬ 
cept we used r = 0.4, a = 0 along the entire interface. 
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TABLE I: Model parameters. We used the same parameters in [27]. Parameters above the horizontal line were also 


used, in the 

same way, in 

[9], 

Name 

Symbol 

Value 

Slider length (x) 

L 

140 mm 

Slider height (z) 

H 

75 mm 

Slider width (y) 

B 

6 mm 

Number of blocks 

W 

57 


w 

31 

Slider mass 

M 

75.6 g 

Block mass 

m 

M/iN^N,) 

Young’s modulus 

E 

3 GPa 

Bulk spring modulus 

k 

WE/4 

Bulk spring length 

1 

L/{N:, - 1) = H/iN, - 1) 

Damping coefficient 

y 

VO.l/cm 

Normal load 

Fn 

1920 N 

Elastic foundation modulus 

kf 

k/2 

Driving spring modulus 

K 

4MN/m 

Driving height 

h 

5 mm 

Driving speed 

V 

0.4 mm/s 

Threshold force coefficient 

f^s — /thres/ /v 0.4 

Slipping force coefficient 

— /slip/ /v 

0.17 

Number of interface springs per block Ns 

100 

Interface spring stiffness 

kij 

y39.2 GN/m2/iv,y 

Slipping time mean 

tR 

2 ms 

Slipping time standard deviation 

Stn 

0.6 ms 

Triggering region width 

^trigger 

22.5 mm 

Triggering region prestress 

^trigger 

0.3 

Time step duration 

At 

2-10"'^s 

Extra damping coefficient 

a 

7?/40 


Appendix D: Slow fronts for systematical variation 
in junction force distributions 

In this appendix we present the data underlying 
Fig. 13c in the main text. 

Figure 21 shows the full evolution of the events corre¬ 
sponding to a range of junction force distribution widths 
(7 and with tq = 0.05. For the lowest a (strong inter¬ 
faces) the front arrests, which is also reflected in Fig. 9. 
For the fronts that do reach the leading edge, the fronts 
are gradually changing from slow to fast propagation. 
These events plot in the “Slow” region of the front type 
phase diagram in Fig. 9 because they would arrest if the 
slow front mechanism was turned off. For each simu¬ 
lation we define the extent of the slow front region by 
visual inspection of the panels in Fig. 21. We then mea¬ 
sure the front speed in the same way as usual, described 
in Appendix A. From these data we determine the aver¬ 
age, minimum and maximum slow front speed values for 
these events. 
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FIG. 21: Color) Underlying data for Fig. 13c: the 
events that lie at r = 0.05 in Fig. 9, shown as in Fig. 3. 

From top left to bottom right, the width of 0 (/t) is 
c^/(/thres ~ /slip) = 0.059, 0.088, and then 0.117 to 0.528 
in steps of 0.059; the event for cr = 0 is not shown, but 
it too arrests early on. For each panel, the extent of the 
slow front was determined by visual inspection. On the 
left, the slow front starts after the step (ca. 

X = 0.04 m). On the right, because the transition to a 
fast front is more gradual, the boundary was set where 
the trend of near constant front velocity is broken. The 
exact limits chosen can be seen from the horizontal 
extent of the data in Fig. 13c. Panels (a-c) show a 
larger time interval than the rest of the panels. No data 
for Fig. 13c were extracted from the bottom right panel: 
even though the front has a slow part and plots in the 
slow part of the front type phase diagram, the slow part 
is too intermixed with the transition to fast rupture for 
the analysis to be carried out. 
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